Effect of Nutrient Enrichment and Turbidity on Interactions Between Microphytobenthos and a Key Bivalve: Implications for Higher Trophic Levels

Benthic diatoms are a high-quality food resource providing essential fatty acids to benthic grazers. Different stressors may alter the proportion of diatoms and other microalgae and thus can affect the quality as well as quantity of food available to benthic consumers. Microphytobenthos (MPB) lipid biomarkers were assessed in a field experiment to elucidate changes to the biosynthesis of fatty acids (FA) under nitrogen (N) enrichment (three levels) at eight intertidal sites that spanned a turbidity gradient. Influences on the flow of carbon and energy were determined using FA biomarkers of a functionally important deposit-feeding tellinid bivalve (Macomona liliana). Site-specific effects of N enrichment were detected in MPB quantity and quality measurements. Enrichment generally increased MPB biomass (chl a) across all sites, while the proportion of diatom associated fatty acid biomarkers was more variable at some sites. Analysis of sediment FA biomarkers and environmental variables suggested that changes to the microbial community composition and quality were related to water clarity and mud content of the bed. The ability of the MPB to utilize the increased nitrogen, as indicated by the resource use efficiency index, was also important. Despite the increase in MPB biomass, lipid reserves in the tissue of M. liliana, a primary consumer of MPB, were reduced (by up to 6 orders of magnitude) in medium and high N addition plots compared to control plots. Further, the nutritional quality of the bivalves to higher trophic levels [indicated by a lower ratio of essential FAs (ω3:ω6)] was reduced in high treatment plots compared to control plots suggesting the bivalves were adversely affected by nutrient enrichment but not due to a reduction in food availability. This study suggests anthropogenic nutrient enrichment and turbidity may indirectly alter the structure and function of the benthic food web, in terms of carbon flow and ecosystem productivity. This may indirectly change the interactions between MPB and key bivalves as suspended sediment concentrations and nutrient enrichment continue to increase globally. This has implications for various ecosystem functions that are mediated by these interactions, such as nutrient cycling as well as primary and secondary production.

Benthic diatoms are a high-quality food resource providing essential fatty acids to benthic grazers. Different stressors may alter the proportion of diatoms and other microalgae and thus can affect the quality as well as quantity of food available to benthic consumers. Microphytobenthos (MPB) lipid biomarkers were assessed in a field experiment to elucidate changes to the biosynthesis of fatty acids (FA) under nitrogen (N) enrichment (three levels) at eight intertidal sites that spanned a turbidity gradient. Influences on the flow of carbon and energy were determined using FA biomarkers of a functionally important deposit-feeding tellinid bivalve (Macomona liliana). Site-specific effects of N enrichment were detected in MPB quantity and quality measurements. Enrichment generally increased MPB biomass (chl a) across all sites, while the proportion of diatom associated fatty acid biomarkers was more variable at some sites. Analysis of sediment FA biomarkers and environmental variables suggested that changes to the microbial community composition and quality were related to water clarity and mud content of the bed. The ability of the MPB to utilize the increased nitrogen, as indicated by the resource use efficiency index, was also important. Despite the increase in MPB biomass, lipid reserves in the tissue of M. liliana, a primary consumer of MPB, were reduced (by up to 6 orders of magnitude) in medium and high N addition plots compared to control plots. Further, the nutritional quality of the bivalves to higher trophic levels [indicated by a lower ratio of essential FAs (ω3:ω6)] was reduced in high treatment plots compared to control plots suggesting the bivalves were adversely affected by nutrient enrichment but not due to a reduction in food availability. This study suggests anthropogenic nutrient enrichment and turbidity may indirectly alter the structure and function of the benthic food web, in terms of carbon flow and ecosystem

INTRODUCTION
Global population increases and land-use change are causing degradation of coastal and estuarine environments (Thrush et al., 2004;Foley et al., 2005). Primarily this degradation is due to the inputs of contaminants, nutrients and sediment from the catchment area that have drastic effects on sediment-dwelling organisms that are essential for ecological function (Thrush et al., 2004). In New Zealand, nearly 200 million tons of terrestrial soil are transported annually from land to sea (Hicks et al., 2011), with sediment deposition smothering, and altering the behavior of functionally important infauna (Cummings et al., 2003;Woodin et al., 2012;Townsend et al., 2014). This in turn can modify biogeochemical gradients and productivity in the sediment (Lohrer et al., 2006;Norkko et al., 2010). Increases in sedimentation are often associated with greater nutrient availability within estuarine systems and reductions in water clarity as turbidity increases (De Jonge et al., 2002;Lovelock et al., 2007). The shift in resource limitation (from nutrients to light) can degrade ecological networks and significantly affect ecosystem function as species interactions are modified (O'Meara et al., 2017).
Primary productivity is the source of energy in most marine food webs. Microphytobenthos (MPB) on the sediment surface can be stimulated by increased nutrient availability associated with the input of fine sediments, but simultaneously become light limited during immersion due to higher turbidity (Drylie et al., 2018). The reduction in light can significantly reduce MPB biomass or photosynthetic capacity (Kromkamp et al., 1995;Cahoon and Safi, 2002;Du et al., 2017), but MPB can also behaviorally and physiologically adapt to low light levels (Cartaxana et al., 2016). Photo-acclimation and migration allow MPB to efficiently adapt to low light conditions but this often results in decreased overall production, with concomitant effects on marine food webs (Consalvey et al., 2004;Migné et al., 2007;Jesus et al., 2009;Serôdio et al., 2012). In addition to the individual cell responses to light, the MPB community composition may be modified as low light or higher nutrient inputs select for different taxa (Gattuso et al., 2006;Hopes and Mock, 2015). For example, the biomass of primary producers may increase with nutrient inputs, but diversity can decrease in response to other limiting resources (Burson et al., 2018) or grazer diversity (Balvanera et al., 2006). MPB quality, and the biosynthesis and transfer of essential fatty acids (FAs) to consumers is rarely considered when assessing the impacts of anthropogenic stressors on soft sediment ecosystem function (Bachok et al., 2006;Bueno-pardo et al., 2018). Yet primary producer quality, indicated by FAs, can often be more important for primary consumers than food quantity (Guo et al., 2016). Light intensity and nutrient availability are considered to be two of the most important factors that determine algae quality and quantity in lakes and streams (Sterner et al., 1997), but the importance of this has rarely been considered in estuarine environments.
Intertidal MPB underpin multiple ecosystem functions and services (Hope et al., 2019). These photosynthetic organisms are the sole source for omega-3 and omega-6 essential fatty acids (EFA) to higher trophic levels and the synthesis of these essential nutrients can be affected by light and nutrient conditions (Hill et al., 2011). MPB are not only important for the flow of energy to higher organisms, but critical for elemental cycling and habitat formation (Sundbäck et al., 2006;Fricke et al., 2017). While bacteria are the main players in the biogeochemical cycling of nitrogen (Cook et al., 2004b), the MPB oxygenate the sediment surface and provide a labile carbon source for N-cycling bacteria, facilitating processes such as denitrification (Tobias et al., 2003;Cnudde et al., 2015). These effects interact through their associations with infauna and bacteria in the sediment Sundbäck et al., 2010). It is therefore important to understand how light and nutrient conditions alter the quantity and quality of MPB not only in their role as a food resource for higher organisms but also for the potential effects on multiple other ecosystem functions.
Several biochemical traits such as stable isotopes and fatty acid biomarkers can be used to determine the nutritional quality of primary producers and consumers (Galloway et al., 2012;Cnudde et al., 2015;Marzetz et al., 2017). In particular, several studies have demonstrated that the synthesis of fatty acid biomarkers can be influenced by changing light, temperature and nutrient conditions (Hill et al., 2011;Jónasdóttir et al., 2019). The resulting changes to the 'quality' of primary producers in terms of FA synthesis can be more important than primary producer quantity for regulating carbon transfer (Müller-Navarra et al., 2000;Marzetz et al., 2017). While some FAs are ubiquitous to all organisms, several taxonomic or functional groups of MPB can be characterized by the presence, ratios or associations between different FAs (Parrish et al., 2000). For example, FAs such as 16:1ω7 (Palmitoleic acid), eicosapentaenoic acid (EPA, 20:5ω3) and docosahexaenoic acid (DHA, 22:6ω3) (Brett and Müller-Navarra, 1997) are mainly produced by MPB such as diatoms, dinoflagellates and cyanobacteria (Volkman et al., 1998). Diatoms are known to be a highly nutritious, high quality, food resource that are rich in these essential fatty acids (EFAs), with specific FAs such as γ-linolenic acid (GLA, 18:3ω6) indicative of cyanobacteria (Yang et al., 2016). Changes to the quality (in terms of EFA synthesis or taxonomical groups) of this food resource will in turn affect the nutritional quality of higher organisms which cannot efficiently synthesize EFAs de novo (Parrish et al., 2000;Bell et al., 2003). Omega-3 FAs (in particular EPA and DHA) are not only vital for the growth and development of herbivores but are critical for reproductive success in many bivalve species (Fearman et al., 2009) and several other physiological processes in organisms including humans (Knauer and Southgate, 1999;Emata et al., 2004;Institute of Medicine [IOM], 2011;Galloway and Winder, 2015;Sprague et al., 2016).
It has been demonstrated that eutrophication in lake systems can reduce the proportion of EPA and DHA transferred to higher organisms from primary producers (Strandberg et al., 2015;Taipale et al., 2016). Furthermore, it is known that cyanobacteria and chlorophytes often dominate the phytoplankton in eutrophic waters (Scheffer and Van Nes, 2007;Berthold et al., 2018) which contain different FA biomarkers. However, there is a lack of similar investigations in marine systems where the focus has been on lipid classes rather than specific FAs (Pinturier-Geiss et al., 2002). Changes to key EFAs and other FAs associated with nutritious diatoms can therefore be used as primary producer and consumer quality indicators (Antonio and Richoux, 2014), and may help to reveal the effects of various anthropogenic stressors such as eutrophication in marine systems. Elevated nitrogen levels can lead to eutrophication with nitrogen toxicity and hypoxia directly and indirectly affecting the growth, survival and reproduction of organisms such as bivalves (Hickey and Martin, 1999;Camargo and Alonso, 2006;Thrush et al., 2017). However, changes to the quantity and quality of their primary food resource, the MPB may also indirectly affect consumer health, influencing activity levels and reproductive success.
Our main objective was to assess the combined effect of light and nutrients on MPB quantity and their nutritional quality as a basal food resource for a functionally important depositfeeding bivalve (M. liliana) within the natural environment. We manipulated sediment porewater nitrogen concentrations across eight estuarine sites in the North Island of New Zealand that encompassed a natural gradient of turbidity. We hypothesized that porewater nitrogen enrichment would alter both the quantity and quality of MPB and the transfer of energy and nutrients from MPB to the deposit-feeder. A secondary hypothesis was that the effects of elevated porewater nitrogen would be site-dependent due to the differences in turbidity across the sites.

Experimental Design
To assess the effects of elevated porewater nitrogen, we conducted a field experiment across eight sites in four estuaries of the North Island of New Zealand that represent a natural turbidity gradient (Figure 1 and Supplementary Table S1). Photosynthetically active radiation (PAR) sensors (Odyssey, Dataflow Systems Pty Ltd., Christchurch, New Zealand) were deployed 10 cm above sediment surface at each intertidal site for approximately 7 months prior to sampling to measure average daily high and low tide PAR (sensors were exposed during low tide). Daytime, mean high-tide (immersed) PAR ranged from < 150 µmol m −2 s −1 to > 900 µmol m −2 s −1 across the selected sites (Figure 2A). As a proxy for turbidity, we derived a percentage clarity measure based on the proportion of light detected from exposed PAR sensors (during daytime low tide) and the submerged sensor (during daytime high tide), as a relative clarity index ( Figure 2B). A clarity value near 100% indicates high water clarity/low turbidity (as clear water allows more of the incident light to the bed during immersion) and a value closer to 0% indicates low water clarity/high turbidity.
At each site, three replicate 9 m 2 plots were randomly allocated to one of three nitrogen treatments; high fertilizer addition (600 g N m −2 ), medium fertilizer addition (150 g N m −2 ), or a disturbance control (0 g N m −2 ). Previous N enrichment experiments have successfully simulated eutrophic estuarine sediments using the above concentrations and application method (see below, Douglas et al., 2016). As nitrogen is typically the limiting nutrient in these systems (Rodil et al., 2011) a slow release nitrogen only fertilizer (Nutricote R N fertilizer [140 days, 40-0-0 N:P:K]) was added to 180 evenly distributed holes (20 holes m −2 ) within each 9 m 2 plot to 15 cm depth using a 3 cm diameter hand-held corer. Douglas et al. (2016) have demonstrated that burial at this depth successfully elevates nitrogen in surface sediments (top 0-2 cm). Nutricote R quickly hydrolyzes to ammonium (NH 4 + ) from urea (Lomstein et al., 1989). Ammonium (NH 4 + ) is the most common form of nitrogen in New Zealand soft sediments (Tay et al., 2013) and the form of N most readily assimilated by microalgae (Underwood and Kromkamp, 1999;Risgaard-Petersen et al., 2004;Riekenberg et al., 2017). Control plots at each site (n = 3) were cored and the sediment then replaced without the addition of fertilizer to replicate the physical disturbance associated with the fertilizer additions in other treatments. All treatments were randomly distributed within the area.

Rationale for Species and Indicators Used in the Study
Macomona liliana is a facultative deposit-feeding tellinid that plays a key role in multiple ecosystem functions on New Zealand intertidal sandflats. Their activities can alter rates of organic matter degradation and nutrient cycling Volkenborn et al., 2012;Woodin et al., 2012). These roles are a result of their complex feedbacks with the MPB which they graze upon and stimulate through nutrient regeneration, and facilitation of biogeochemical cycling through their sediment reworking. M. liliana were used as a model benthic primary consumer while intertidal sediment surface (0-2 cm) samples were used to characterize MPB quality and quantity. While M. liliana are capable of suspension feeding, their primary feeding mode is deposit-feeding (Pridmore et al., 1991) and strong but complex relationships with the MPB have been demonstrated (Thrush et al., 2004). Furthermore, previous investigations have demonstrated that resuspended MPB are the primary food resource in NZ estuaries for both deposit and suspension feeders (Savage et al., 2012;Jones et al., 2017). As FAs are a primary component of lipids and play a  Supplementary Table 1. fundamental role as an energy and nutrient source to higher organisms, FA biomarkers and other biochemical traits (lipid content, MPB biomass) were employed to determine MPB quantity and quality, MPB community responses and dietary intake in M. liliana.
As individual fatty acid biomarkers provide limited information, key ratios indicative of MPB community composition, were selected for examination ( Table 1). These included: the diatom index of the MPB assemblage to determine the relative dominance of diatoms over dinoflagellates (Antonio and Richoux, 2014); and the proportion of bacterial FAs (BaFA) which can indicate the proportion of aerobic and anaerobic bacteria and degree of decomposition (Jaschinski et al., 2008;Ruano et al., 2012). These indices were also used to determine the relative ingestion or retention of essential FAs in bivalves or the ingestion of carbon that has first been processed through the microbial loop. Other FA biomarkers such as the ratio of ω3:ω6 PUFAs, and the proportion of EPA + DHA, were also examined as useful indicators of the nutritional quality of the primary producers, or 'health' of the bivalves (Vargas et al., 1998;Ruano et al., 2012).

Sample Collection and Processing
Sampling was conducted between October and November 2017, approx. 7 months after the initial nitrogen enrichment of the plots. During mid-morning low tides, four sediment cores (0-1 cm depth, 2.6 cm dia.) were randomly collected and pooled for biochemical analyses from each plot (giving 3 replicates per nutrient treatment, per site). An additional 5 cores were sampled and pooled from each replicate plot to determine particle size distribution and porosity (0-2 cm depth, 2.6 cm dia.) and 4 cores were pooled from each replicate plot for porewater nutrient analysis (0-2 cm depth, 2.6 dia.); again providing 3 replicates per nutrient treatment per site. The pooled sediment cores from within each replicate plot were used to ensure adequate sediment for subsequent analyses, and to account for the spatial variability across the 9 m 2 plots in biochemical properties, such as microphytobenthic (MPB) biomass.

FIGURE 2 | (A)
Mean daily photosynthetically active radiation (PAR, ± SE) measured 10 cm above the sediment surface for daytime low tide (white bars) and high tide (gray bars) for the 7 months period between the experiment setup and sampling. Daytime averaged from 2 h after sunrise to 2 h before sunset on days when low tides corresponded with daylight hours. (B) Relative clarity at each site over the 7 month period. Values close to 100% indicate the majority of incident light (emerged sensors) is detected by the submerged sensor (during immersion periods). Values close to zero mean very little incident light is reaching the submerged sensor.
Sediment samples for biochemical analysis (see below) were foil wrapped, immediately flash frozen in liquid nitrogen and transported to a −80 • C freezer until processed to preserve FA biomarkers and chlorophyll content. Porewater samples were stored on ice in the dark until extracted in the laboratory. Duplicate 13 cm dia. cores were collected to a depth of 15 cm for the gathering of adult M. liliana specimens (> 20 mm shell length) from within each plot. Specimens were rinsed of sediment, foil wrapped and flash frozen in the field to be transferred to a −80 • C freezer until further analysis. M. liliana were collected from all sites with the exception of MH2 and WR1, where specimens could not be detected in plots across all three nutrient treatments. These sites were excluded from further FA analysis of M. liliana tissue. Particle size distribution (PSD) was determined from freeze-dried sediments using a Malvern Mastersizer 2000 (range 0.05-2000 µm) (Singer et al., 1988) after the removal of organic matter with 10% hydrogen peroxide solution and disaggregation using 6% calgon.

Biochemical Properties
Total M. liliana biomass (incl. shell) was measured from the cores in each plot and biomass standardized to kg M. liliana m −2 . Whole organisms (n > 3) were freeze dried and homogenized to determine the animals' lipid reserves (conditional status) and fatty acid biomarkers. Sediment OM content was determined by loss on ignition of dried sediments at 450 • C for 4 h, following Parker (1983). Porewater was extracted by adding 4 mL of deionized water to sediments, centrifuging at 3500 RPM for 10 min and filtering the extract through Whatman GF/C (1.1 µm) filters (Lohrer et al., 2010). Thawed extracts Indicators of MPB and bivalve nutritional quality are centered around the importance of omega-3 fatty acids (ω3).
were analyzed for ammonium (NH 4 + ), using standard methods for seawater analysis on a Lachat QuickChem 8500 + FIA (Zellweger Analytics Inc., Milwaukee, WI, United States) with concentrations corrected for sediment porosity and dilutions. Porewater concentrations were used to determine the effectiveness of nitrogen enrichment of the plots across the sites. Further, to assess the effects of N enrichment across the turbidity gradient the ratio between each sites high-tide light level (PAR) and the mean porewater NH 4 + (light: NH 4 + ) of treatment plots were calculated. This ratio has previously been used to assess the changes in microalgae fatty acid synthesis due to shifts in the ratio of available light and nutrients (Hill et al., 2011).
Chlorophyll a (chl a) and pheophytin (pheo) contents were determined from surficial sediments using a 90% acetone extraction to quantify MPB biomass and degradation products respectively (Lorenzen, 1967). Pigments were extracted in the dark for 24 h at 4 • C using homogenized, freeze-dried sediment (∼1 g). Acetone was selected as it is best suited for spectrophotometric assay (Ritchie, 2006). Changes in resource use efficiency (RUE) was calculated as the ratio between sediment chlorophyll a concentrations and porewater concentration of the limiting nutrient; nitrogen in these systems. RUE is an ecological measure that defines the proportion of supplied resources that is then converted into new biomass (Hodapp et al., 2019). It is a useful metric to illustrate the change in resource use efficiency as higher concentrations of key nutrients become available to the MPB (Eriksson et al., 2017). An increase in RUE indicates that the MPB are efficiently assimilating the identified nutrient (nitrogen in this case) and converting it to biomass. Aliquots of tissue from freeze-dried bivalves and sediments were used to assess total lipids using a modified Bligh and Dyer (1959) and lipid content determined using the sulfo-phospho-vanillin (SPV) spectrophotometric method (Byreddy et al., 2016). Total fatty acid (TFA) composition was determined following a one-step direct transesterification method (Lepage and Roy, 1986;Zárate et al., 2016). The full details of the extraction and identification of FAs is provided in the Supplementary Material. All biochemical properties discussed in this section were standardized by either sediment/tissue weight or by the surface area of the sediment core as deemed appropriate (Tolhurst et al., 2005).
Identified FAs were first expressed as a percentage of the total FAs (TFAs) identified in each sample (% TFA) or as a ratio between two FAs if a dominance indicator was used. All FAs are designated as X:YωZ, where X is the number of carbon atoms, Y is the number of double bonds and Z is the position of the ultimate double bond from the terminal methyl.

Data Analyses
Changes in MPB biomass and key FA biomarkers were tested across the different nutrient treatments (control, medium, high). To test the effects of nutrient enrichment and site on biochemical properties and FA biomarkers of the sediment and bivalves ( Table 2), separate two-way PERMANOVAs based on Euclidean distance matrices (PRIMER, V7, Anderson et al., 2008) were run using a fixed factor nutrient (3 levels), a random factor of site (8 levels) and an interaction term. A total of 9999 permutations were applied to residuals under a reduced model. Monte Carlo p-values were also considered for the significance of differences among factors. Pairwise comparisons were performed using a post hoc test with the Bonferroni correction.
To summarize relationships between MPB quality indicators and sediment properties ( Table 2) we used canonical analysis of principal components (CAP, Anderson et al., 2008), based on the Euclidean distance matrix of MPB quality indicators. Similarly, relationships between M. liliana biomarkers and biochemical traits (listed in Table 2), environmental and MPB/sediment indicators were summarized by CAP analysis based on the Euclidean distance matrix of M. liliana FA biomarkers. CAP was used rather than distance-based redundancy analysis (dbRDA) as the feedbacks between the different variables precludes the definition of any clear explanatory and response variables. All data used in the CAP analyses were normalized using a fourthroot transformation. The number of PCO axes (m) was chosen to exclude redundant axes, but include as much of the original variability in the dataset while minimizing the leave one out residual sum of squares (Anderson et al., 2008). X Quality indicators marked with an X were used in multivariate analyses. Environmental properties were assessed against these matrices for both sediment and bivalve data. Additional indicators from the sediment were explored against the bivalve data to assess whether changes in the quality and quantity of MPB in the sediment were reflected in changes to the quality indicators of the bivalves. Diatom index under quality indicators is the diatom index of bivalves.

Sediment Nutrients and Grain Size
The relative clarity of overlying water varied from 19% in the Manukau Harbor (MK, low clarity and high turbidity) to 93% at Whangarei 2 (WR2, high clarity and low turbidity, Figure 2B). The addition of fertilizer to the sediment resulted in elevated surface porewater NH 4 + concentrations (0-2 cm sediment depth) across all sites. Enrichment was proportional to the amount of fertilizer added, in that the final elevated porewater NH 4 + concentrations varied between sites ( Figure 3A), but higher N additions elevated the final concentrations of medium and high plots compared to controls (C < M < H plots, Pseudo-F = 9.0, P < 0.001). Mean grain size ranged from very fine to medium sand (63-500 µm; Wentworth, 1922), and all sites had a mean mud content < 10% (Figure 3B). Neither mean grain size nor mud content, varied significantly with the nutrient enrichment of the sediment, but site differences were detected for both factors (Pseudo-F = 25.7, P < 0.001, Pseudo-F = 20.7, P < 0.001 respectively, Table 3).

Changes to Microbial Communities
Nitrogen enrichment and site differences in turbidity influenced various microbial indicators in the sediment ( Table 3). Nitrogen enrichment typically increased MPB biomass (chl a) but effects were site-specific (Figure 4A). At a number of sites, [MH1, MH2, WR1, WR2], high variability in the MPB biomass prohibited the detection of significant increases between control plots and medium or high N addition plots, while at two of the sites [MH1, WR1], MPB biomass was elevated in medium treatment but not in high treatment plots ( Figure 4A). Several sites also exhibited a higher sediment diatom index [MK, WR3, MH3, MH1 and control-medium treatments at WR1], while the other sites displayed no differences between N treatments ( Figure 4B). Interestingly, the proportion of diatoms in the MPB (diatom index) of control sediments across the sites exhibited no significant relationship with water clarity (P < 0.05). The EPA content of the sediment (characteristic FA biomarker for diatoms), increased with nitrogen treatment but there were no significant site differences detected ( Table 3). The MPB were not effectively utilizing the higher N concentrations available in the porewater in enrichment plots as indicated by the reduced resource use efficiency (RUE, Figure 5). The proportion of aerobic and anaerobic bacteria-associated FAs (BaFAs) in the sediment varied primarily across the sites and not specifically with elevated nitrogen.
We further examined Pearson correlation coefficients between the turbidity proxy (% clarity) and measured variables for each nutrient treatment independently. With the exclusion of the Manukau (MK) site, MPB biomass of control plots was higher at sites with greater water clarity (r 2 = 0.73, Figure 6). The relationship between MPB biomass and clarity essentially disappeared with N enrichment (r 2 = 0.08 and 0.07 for medium and high N treatments respectively). Similarly, while BaFAs in control plots correlated significantly with water clarity (Figure 7, r 2 = 0.67), the relationships in medium (r 2 = 0.05) and high (r 2 = 0.05) nitrogen treated plots were diminished.

Changes to the Bivalve M. liliana
The biomass (Figure 8) and total lipid content ( Figure 9A) of M. liliana was reduced with nitrogen enrichment at the majority of sites, with a significant interaction between site and N enrichment detected (Table 3). While most sites had significant reductions in both metrics, no effect of nitrogen enrichment was apparent for biomass at MK (Figure 8) or lipids at WG ( Figure 9A). Coinciding with a decrease in the overall biomass of bivalves in N enriched plots and a depletion in lipid reserves of the bivalves that were collected in the plots, there was a significant decrease in the nutritional quality of the bivalves (ratio of ω3:ω6 fatty acids (Figure 9B) across the majority of sites with nitrogen additions ( Table 3). The effects of elevated nitrogen on the diatom index and the EPA + DHA content of M. liliana tissue was also site dependent. At Manukau (MK) and Whangarei site 3 (WR3) the proportion of diatoms and intake of EPA + DHA in the diet of M. liliana were significantly reduced under elevated nitrogen levels (Figures 9C,D) whereas at Whangateau Harbor (WG), the intake of diatoms (diatom index), bivalve quality (ω3:ω6 ratio and the EPA + DHA content) all increased (Figures 9B-D).

Relationships Between M. liliana, MPB and Environmental Variables
Overall the CAP trace statistic (tr = 0.68, P < 0.001) together with the ordination plot suggests that the matrix of MPB quantity (biomass) and quality (FA biomarkers) in the sediment were strongly related to several measured sediment and environmental Mean mud content (%, ± SE, n = 3) for each treatment at each site, and nutrient treatment. Nutrient addition treatments; white bars -control plots (0 g N m -2 addition), light gray -medium plots (150 g N m -2 addition and dark gray -high plots (600 g N m -2 addition).
properties (Figure 10). Diagnostics revealed the first three axes encapsulated around 96% of the variability in the resemblance matrix and minimized the residual sum of squares. The variables most strongly related to CAP1 were the sediment mud content (r 2 = 0.47), the resource use efficiency (RUE) of the MPB (r 2 = 0.25) and the percent clarity (r 2 = -0.17). The second axis (CAP2) was related to porewater NH 4 + concentration (r 2 = 0.36) and the ratio between light and nutrient availability (Light: NH 4 + , r = 0.20). The third axis contributed even less (CAP3, δ 2 = 0.11) and was primarily related to the Light: NH 4 + ratio and RUE of the MPB (r 2 = 0.15 for both).
Interestingly, the CAP ordination of M. liliana tissue data displayed contrasting results (Figure 11, tr = 0.76, P < 0.01). Instead of being dominated by one axis, the first two axes were of similar importance (δ 2 = 0.32 and 0.30, respectively), with the third axis contributing markedly less (δ 2 = 0.09). The environmental and MPB variables associated with the first axis (CAP1) were the mean grain size of the bed (r 2 = −0.45), however the 2nd axis (CAP2) was associated with the ratio between the light availability and porewater nitrogen (Light:NH 4 + ratio, r 2 = −0.35) and the proportion of diatoms (r 2 = 0.38) and bacterial-associated FAs (r 2 = −0.22) as well as the clarity of the water column (r 2 = −0.14, Figure 11).

DISCUSSION
This study demonstrates that elevated porewater nitrogen affects the proportion of diatom-and bacteria-specific fatty acid biomarkers available in the sediment to primary deposit feeding consumers. Changes to FA biomarkers and MPB biomass were different across the sites, and a number of variables were related to the overlying water clarity. Our data suggest that the proportion of bacterial FA biomarkers was higher when water column turbidity was lower. Further, the positive correlation between MPB biomass and water clarity under ambient nitrogen conditions (controls) suggests reduced light availability in turbid estuaries may limit MPB growth. The exception to this was of course the high biomass MPB detected at the Manukau harbor (MK) despite high turbidity. Low light and high nutrient availability has been evident at the MK site for a number of years (Cahoon and Safi, 2002) and as such the MPB present may be photoacclimated to lower light (Veach and Griffiths, 2018).
Previous studies have demonstrated that productivity during low tide can compensate for the lack of productivity during immersion in highly turbid sites similar to the Manukau (Drylie et al., 2018). It has previously been highlighted that there is a need for multiple stressor studies that incorporate a wide range of habitat conditions and species interactions (Crain et al., 2008). Our context-dependent results emphasize the complexity of these interactions and the fundamental need for large-scale field studies to encapsulate a wide range of ecological effects. Our findings, suggest that effects of increased nitrogen loading on the functional role of MPB may depend on the environmental context (background nutrient levels and sediment type). Context dependency and legacy effects in ecological field studies are well documented, but can provide a wealth of information that cannot be provided by laboratory investigations (Thrush et al., 2000;Norkko et al., 2006). In these predominantly noncohesive, low nutrient soft-sediments, the addition of nitrogen may have increased the quantity and quality of the MPB but their functional role in the system and trophic interactions were altered as the increase in essential FAs associated with diatoms were not reflected in the consumers The quality and quantity of MPB both increased with nitrogen enrichment but the magnitude of the effects were again site dependent. The correlations between MPB biomass and water clarity observed in control plots disappeared with N enrichment of the sediment suggesting the increase in nitrogen availability in the sediment facilitated the growth of the MPB by providing them with the nutrients required during low tide periods. Primary producers of course require both light and nutrients for growth and to reproduce (Kromkamp et al., 1995). Carbon is obtained through photosynthesis which relies on light availability, while key nutrients such as nitrogen and phosphorus are obtained from their environment (Longphuirt et al., 2009). The positive response of the MPB biomass to nitrogen additions suggests that the MPB growth across the sites before experimental enrichment was nitrogen limited (Underwood and Kromkamp, 1999). Non-cohesive sediments, associated with greater hydrodynamic energy, often have lower MPB biomass and nutrient capacity than muddy sediments (Cook et al., 2004a;Eyre et al., 2011). Remineralization and tight coupling of the nutrients within sandy sediments can however provide sufficient nutrients to sustain the MPB community and benthic metabolism (Cook et al., 2004a;Sundback and Mcglathery, 2005). Nonetheless, these differences in the sedimentary environment influence the rate and pathway of OM degradation (Ferguson et al., 2004). While an increase in RUE would suggest that the MPB are effectively using the additional nitrogen added to the system, a reduction (which we observed) suggests that N is now in excess following porewater enrichment and the MPB demand is no longer able to keep up with the supply of this nutrient from the porewater.
As inputs of nitrogen and sediment from catchments to coastal and estuarine environments continues to increase (Howarth, 2008;Hicks et al., 2011), these two stressors synergistically interact to affect the quality and quantity of primary food resources and their consumers. Elevated sediments enhance turbidity and limit light availability to the MPB and decrease their ability to assimilate nutrients at the sediment-water interface (Underwood and Kromkamp, 1999;Pratt et al., 2014), as we have seen in the RUE above. This causes an increase in the release of nutrients to the water column that further stimulate the phytoplankton enhancing turbidity effects. While we have examined the effects of N enrichment across a natural turbidity gradient that is attributed more to suspended sediment concentrations rather than phytoplankton blooms, our data suggest that as nutrient concentrations increase there may be a breakdown in the relationship between key consumers such as M. liliana and MPB across the light gradient. Understanding and predicting the effects of multiple stressors in the field is difficult, but critical to improving our understanding of anthropogenic effects in the real world (Crain et al., 2008;Thrush et al., 2012). Our study provides real-world information on the complexity associated with trophic interactions across a range of different sedimentary habitats. Furthermore, we provide further evidence of the potential impacts of eutrophication in understudied, low nutrient systems on carbon flow from primary producers to consumers. We stress that this is not only dependent on the adaptation of the MPB to particular stressors but will be influenced by feedbacks associated with differences in the response of higher trophic levels to the same perturbations Crain et al., 2008).
At two sites (MH2 and WR1), M. liliana were absent in nitrogen enriched plots. At sites where we did obtain M. liliana specimens, these bivalves from enriched plots had decreased lipid energy reserves. Plausible explanations for the reduced lipids energy reserves and lower nutritional quality of the bivalves include the accumulation of ammonium in the sediment becoming a chemical stressor to biota (Camargo and Alonso, 2006;Simpson et al., 2013). In other bivalves this type of stress response, the reductions in lipid content, has been related to oxidative stress and lipid peroxidation (Coutteau et al., 1996;De Almeida et al., 2007).

FIGURE 4 | (A)
Mean chlorophyll a concentration and (B) mean diatom index of the sediment surface, across all sites and nutrient treatments (± SE, n = 3). Nutrient treatments; white bars -control plots (0 g N m -2 addition), light gray -medium plots (150 g N m -2 addition), dark gray -high plots (600 g N m -2 addition).
The diatom-associated FA biomarkers and overall lipid energy reserves of M. liliana were reduced by N enrichment suggesting that grazing on the MPB was diminished in our experimental plots. These changes are consistent with changes in dietary intake, as more diatom rich FAs were available in the sediment (Bell et al., 2003). This reduction in FA intake was not observed in M. liliana inhabiting Whangateau estuary. Here, M. liliana appeared to benefit and increase their uptake of nitrogenstimulated MPB (diatom associated FAs). This is likely due to the relatively shallow, ebb-dominated nature of the estuary (Grace, 1972) and the associated advection of excess organic matter offshore (Middelburg and Herman, 2007). The tight recycling of nutrients in the bed allowed this low-turbidity system to remain productive and healthy. The depletion of lipid reserves can also result in a shift between different biomarkers, as specific FAs are differentially metabolized (Budge et al., 2001;Dridi et al., 2017). Despite the preservation of the DHA/EPA ratio, the ω3:ω6 PUFA ratios of bivalves were generally reduced in enriched plots. This suggests that bivalve nutritional quality was compromised (Jaschinski et al., 2008;Ruano et al., 2012). The higher proportion of bacterial FAs in M. liliana tissue suggests M. liliana ingested more bacteria or organic matter that had been processed via the microbial loop (Mayzaud et al., 1990;Gonçalves et al., 2017). The higher proportion of diatoms in the sediment and the shift in bivalve quality suggests the changes to the MPB community stimulated bacterial communities, which in turn lowered the quality of the bivalves. The shift in ω3:ω6 PUFA ratios is significant, as the transfer of ω3 PUFA, in particular EPA and DHA, is essential for the growth and reproduction in higher trophic levels and supports the provision of food and goods (fish, shellfish and Omega-3) from our coastal marine ecosystems (Kharlamenko et al., 2008;Twining et al., 2016).
Changes in the quality and quantity of food resources to secondary consumers has implications for trophic food webs, as intertidal bivalves are the major prey of many crustaceans, fish, and wading birds. As the condition of these bivalves FIGURE 5 | Mean resource use efficiency (RUE), i.e., the ratio between chlorophyll a concentration and the sediment porewater concentration of the limiting nutrient [N in these systems (RUE calculations based on Eriksson et al., 2017)]. Nutrient treatments; white bars -control plots (0 g N m -2 addition), light gray -medium plots (150 g N m -2 addition), dark gray -high plots (600 g N m -2 addition).
is altered their functional roles as prey for higher trophic levels and as grazers and stimulators of the MPB will also be altered . While nutrient regeneration mediated by bioturbators such as M. liliana may become less influential on the MPB under higher external nitrogen loads, the utilization of nutrients by the MPB community would also likely decrease. This highlights the potential for top down processes to influence the MPB under N enrichment: Diatoms are proliferating under higher nitrogen concentrations not only due to the increase in nitrogen required for growth and reproduction but due to a decrease in grazing pressure. As FAs are essentially incorporated unaltered into the lipids of first order consumers, their FA profiles should reflect the FA profiles of the food they consume (Kharlamenko et al., 2001;FIGURE 7 | Relationship between the water clarity and the mean proportion of bacteria-associated FAs (C15:0 + C17:0) in the sediment across ambient nitrogen conditions (control plots) in all sites. (r 2 = 0.67).
FIGURE 8 | Mean biomass (incl shell) of M. liliana specimens (20 -30 mm length, kg dry weight bivalves m 2 , ± SE, n = 3) as a function of site and nutrient enrichment treatment. Nutrient treatments; white bars -control plots (0 g N m -2 addition), light gray -medium plots (150 g N m -2 addition), dark gray -high plots (600 g N m -2 addition). Dalsgaard et al., 2003;Lebreton et al., 2011). Our data suggests the proportion of diatom-associated FAs in the bivalves decreased with nitrogen enrichment in the majority of sites alongside a concomitant increase in bacterial associated FAs. This could be due to the bivalves feeding less on the MPB due to N enrichment of the sites, or an increase in the metabolism of diatomassociated FAs. Direct chemical stress from other pollutants such as pesticides, pH, temperature and sedimentation have been observed to decrease feeding rates and nutrient uptake in several consumers (Patil, 2011;Wang et al., 2015;McCartain et al., 2017), which would in time reduce their biomass as fat stores are used. As bivalves tend to increase their metabolic activity to survive extended periods of stress from pollutants (Smolders et al., 2004;De Almeida et al., 2007;Patil, 2011) this would further reduce biomass and lipid energy and FA reserves and potentially negatively impact their reproductive output.
Multivariate analysis of relationships between sediment FA biomarkers and environmental variables suggested that changes to the microbial community composition and quality was related to site-specific conditions such as the water clarity and subtle changes in the mud content of the sandy sediments. The effects of nutrient enrichment on FA biomarkers was modulated differently across the sites depending on the light availability. The functional role of intertidal microphytobenthos (MPB) in estuarine food webs have traditionally been characterized by their biomass, contribution to primary productivity and their influence on bacteria and nitrogen pathways at the sediment-water interface (MacIntyre et al., 1996;Serôdio and Catarino, 1999;O'Meara et al., 2017;Koedooder et al., 2019;Vaz et al., 2019). Our findings support previous studies where shifts in the MPB community were apparent with nutrient additions (Pinckney, 1995;Piehler et al., 2010). However, our experiment also reveals that the biosynthesis and transfer of essential FAs to consumers are rare but important components of the impacts of anthropogenic stressors on ecosystem function (see also Bachok et al., 2006;Bueno-pardo et al., 2018). Higher nitrogen inputs into low nutrient systems such as the ones we have studied herein can therefore alter the functional role of MPB as primary producers and as a basal food resource.
Across different sedimentary habitats resuspended MPB contribute substantially to water column productivity (Underwood and Kromkamp, 1999;Jones et al., 2017) and can even comprise up to 70% of the diet of harvested and farmed suspension feeders such as oysters, mussels and cockles (Sauriau and Kang, 2000;Dubois et al., 2007;Morioka et al., 2017). The MPB therefore support various coastal fisheries (Kritzer et al., 2016) and changes to the biosynthesis and transfer of essential FAs (EFAs) to key consumers may therefore have profound impacts on ecosystem functions through altering the functional role of the MPB and the biodiversity of fauna depending on their tolerance to high nutrient concentrations. Our investigation suggests that a key deep-dwelling tellinid bivalve that dominates FIGURE 10 | Canonical analysis of principal components ordination of the quantity and quality indicators of the MPB and their relationship to environmental variables. Canonical axis 1 (δ 2 = 0.41) and canonical axis 2 (δ 2 = 0.15). PW NH4 + -porewater NH 4 + concentration (0 -2 cm depth). MC -mud content of the bed. RUE -resource use efficiency of the MPB (chl a:NH 4 + ratio).Light:NH 4 + ratio -The ratio between high-tide light and porewater NH 4 + concentration.% clarity -The proportion of low tide light that reaches the PAR sensor during high tide, used as a proxy of turbidity. NZ estuaries, is negatively affected by nitrogen enrichment, and this could alter their role in soft sediment functions and their interactions with MPB and bacteria.

DATA AVAILABILITY STATEMENT
The datasets generated for this study are available on request to the corresponding author.

AUTHOR CONTRIBUTIONS
The New Zealand National Science Challenge Sustainable Seas, Dynamic Seas, Tipping Points project (CO1 × 1515) was responsible for the design and implementation of the larger field experiment, this involved many researchers in the implementation, maintenance and sampling of the experiment (see below). JAH conceived the manuscript and produced the first draft of the manuscript. JAH collected and processed fatty acid samples. JAH analyzed data with advice from JH. All authors contributed to the ideas presented in this manuscript, drafts of the manuscript and gave final approval for publication.

ACKNOWLEDGMENTS
The field campaign was conducted as part of a larger 'Tipping points' project funded through the National Science Challenge Sustainable Seas, which involved all authors and a nationwide group of collaborators to which we gratefully thank for the conception of the field experiment, the set-up and maintenance of the nitrogen enriched plots and deployed sensors throughout New Zealand. We would also like to thank the Seafloor Ecology Group at the University of Auckland for sample collection from the Auckland sites, Steph Mangan for the summary of PAR data and NIWA, Hamilton for processing porewater nutrient samples. We wish to thank Erica Zarate and Saras Green of the Mass Spectrometry Centre, Auckland Science Analytical Services. The University of Auckland, New Zealand for assistance with analysis of FAMEs by GC-MS, and to Raphael Bang for data processing.