Assessing biomass and primary production of microphytobenthos in depositional coastal systems using spectral information

In depositional intertidal coastal systems, primary production is dominated by benthic microalgae (microphytobenthos) inhabiting the mudflats. This benthic productivity is supporting secondary production and supplying important services to humans including food provisioning. Increased frequencies of extreme events in weather (such as heatwaves, storm surges and cloudbursts) are expected to strongly impact the spatiotemporal dynamics of the microphytobenthos and subsequently their contribution to coastal food webs. Within north-western Europe, the years 2018 and 2019 were characterized by record-breaking summer temperatures and accompanying droughts. Field-calibrated satellite data (Sentinel 2) were used to quantify the seasonal dynamics of microphytobenthos biomass and production at an unprecedented spatial and temporal resolution during these years. We demonstrate that the Normalized Difference Vegetation Index (NDVI) should be used with caution in depositional coastal intertidal systems, because it may reflect import of remains of allochthonous pelagic productivity rather than local benthic biomass. We show that the reduction in summer biomass of the benthic microalgae cannot be explained by grazing but was most probably due to the high temperatures. The fivefold increase in salinity from January to September 2018, resulting from reduced river run-off during this exceptionally dry year, cannot have been without consequences for the vitality of the microphytobenthos community and its resistance to wind stress and cloud bursts. Comparison to historical information revealed that primary productivity of microphytobenthos may vary at least fivefold due to variations in environmental conditions. Therefore, ongoing changes in environmental conditions and especially extreme events because of climate change will not only lead to changes in spatiotemporal patterns of benthic primary production but also to changes in biodiversity of life under water and ecosystem services including food supply. Satellite MPB data allows for adequate choices in selecting coastal biodiversity conservation and coastal food supply.


Abstract:
In depositional intertidal coastal systems, primary production is dominated by benthic microalgae (microphytobenthos) inhabiting the mudflats. This benthic productivity is supporting secondary production and supplying important services to humans including food provisioning. Increased frequencies of extreme events in weather (such as heatwaves, storm surges and cloudbursts) are expected to strongly impact the spatiotemporal dynamics of the microphytobenthos and subsequently their contribution to coastal food webs. Within north-western Europe, the years 2018 and 2019 were characterized by record-breaking summer temperatures and accompanying droughts. Field-calibrated satellite data (Sentinel 2) were used to quantify the seasonal dynamics of microphytobenthos biomass and production at an unprecedented spatial and temporal resolution during these years. We demonstrate that the Normalized Difference Vegetation Index (NDVI) should be used with caution in depositional coastal intertidal systems, because it may reflect import of remains of allochthonous pelagic productivity rather than local benthic biomass. We show that the reduction in summer biomass of the benthic microalgae cannot be explained by grazing but was most probably due to the high temperatures. The fivefold increase in salinity from January to September 2018, resulting from reduced river run-off during this exceptionally dry year, cannot have been without consequences for the vitality of the microphytobenthos community and its resistance to wind stress and cloud bursts. Comparison to historical information revealed that primary productivity of microphytobenthos may vary at least fivefold due to variations in environmental conditions. Therefore, ongoing changes in environmental conditions and especially extreme events because of climate change will not only lead to changes in spatiotemporal patterns of benthic primary production but also to changes in biodiversity of life under water and ecosystem services including food supply. Satellite MPB data allows for adequate choices in selecting coastal biodiversity conservation and coastal food supply.   Yes -all data are fully available without restriction

106
The aim of this study is to test generic methods to determine biomass and productivity of microphytobenthos by

123
The study area is the Dollard, the innermost part of the Ems estuary, which is enclosed between the Netherlands

124
in the west and Germany in the east (Fig 1). The Dollard has a surface area of 103 km 2 of which 81% consists of

141
The subsequent increase in turbidity as a result of the loss of these sediment sinks [60] has caused concern and 142 large-scale measures are being put into place including the construction of artificial saltmarshes at the expense Assessing MPB biomass using spectral information of tidal flat systems [61]. Although it appears that the increased turbidity has resulted in a reduction of pelagic 144 primary production [58], the historical and future effects on the benthic production remains unclear.

153
The analysis followed a 'biological approach', meaning that organic matter and calcium carbonate were not 154 removed [62]. Mud is defined as the fractions <63 µm (volume), and mud percentages were calculated as the 155 contribution of this fraction to the total volume.

161
Before analysis, the samples were defrosted and 20 ml of 90% acetone was added to them, suspensions were 162 mixed and stored overnight at 4°C in the dark. The next day, samples were homogenized, and 8 ml of sample 163 extract was centrifuged at 3000 RCF for 10 min. From the sample extract, 3 ml was pipetted in a cuvette and 164 measured using a fluorescence spectrophotometer (F-2500 Hitachi). When concentrations were too high to be 165 measured, samples were diluted with 90% acetone. Samples were measured again after adding 2 drops of 10%

167
Following [63], concentrations of uncorrected (not acidified) chlorophyll-a concentrations (CHLa_u; mg m -2 ), 168 corrected chlorophyll-a concentrations (CHLa_c; mg m -2 ) and pheophytin-a (PHEOa; mg m -2 ) were determined 169 using the following equations: In these equations, Rb (unitless) is the fluorescence signal of the sample before adding the acid solution; Ra 174 (unitless) is the fluorescence signal after adding the acid solution, Vextration is the volume acetone added to the 175 sample (ml), Vsample is the sampled area (cm 2 ) and DF is the dilution factor (unitless).

176
The other parameter values (Fs and r) are derived from calibrations that were performed once for each sampling   chlorophyll-a concentration of the slurry (µg l -1 ) was determined as described above (2.2) using triplicates using 191 10 ml of slurry instead of the sediment cores.

192
The incubation procedure followed that of [64], with small modifications. In the radioisotope laboratory, each 193 slurry was well agitated while 2 ml subsamples were pipetted into 20 ml glass incubation vials.

212
For each sample, the corrected rate of disintegrations per minute (dpm) was calculated as the measured dpm of 213 that sample (dpmsample) minus the average dpm of the two dark flasks (dpmdark). DIC is the concentration of 214 dissolved inorganic carbon (mg C l -1 ).

215
To correct for a temperature difference between the in-situ temperature and the temperature during the 216 incubation, a correction factor was applied: The constant 1.05 in equation (2) is a factor to correct for the lower uptake rate of 14 C compared to 12 C, c is a 219 constant with value 1000 to convert units, dpmadded is the dpm as measured in the control bottles corrected for 220 the volume used and t is the duration of the incubation (in hours). The fixation rates were normalised to chlorophyll-a concentrations (from the slurry, see previous paragraph) and with these rates P-E curves were 222 fitted [68-69]. P-E curves were fitted using a model described by [70].

223
To calculate the production of each sampling day, information on irradiance during the sampling day, the light

256
Downwelling irradiance (Es; units W m -2 ) and radiance (Lu; units W m -2 sr -1 ) sensors, respectively were installed 257 on a portable frame and controlled with a field laptop. Spectral data were interpolated at 1 nm (Fig 2). The

258
surface information, contained in Lu, was normalized by Es to derive the surface reflectance: Here, is unitless due to the scaling factor (sr), that accounts for the conversion between radiance and 261 irradiance.

262
The normalized difference vegetation index (NDVI; unitless) was calculated according to: Here  Fig 2). CHLa (+PHEOa) was regressed to NDVI from each of these band settings to obtain sensor-274 specific algorithms.

316
The updated bathymetry of the mudflats, provided by Rijkswaterstaat at 30 m resolution, was resampled to the 317 satellite resolution and used to discriminate the satellite data based on it. Here, we defined two intertidal zones,

318
namely "high" mudflats as those 40 cm or more above NAP and "low" mudflats being those below NAP +40   (Table 2).

336
For each of the sampling periods, the median grain size (MGS) and fraction mud of the sediment were highly 337 negatively correlated (-0.98  r  -0.91; Fig 3).    (Table 3).

391
For all three sampling periods, the sum of corrected chlorophyll-a and pheophytin-a (mg CHLPH_c m -2 ) at

393
for September 2018, also a positive relationship of this sum with pheophytin-a was found (Fig 3).

406
The light-production curves that were fitted for each of the three sampling dates and two stations showed that 407 the slope of the light-limited part of the curve (α β ; mg C (mg CHLa_c) -1 h -1 (PAR µE m -2 s -1 ) -1 ) and the 408 maximum photosynthetic production rate (P β max; mg C (mg CHLa_c) -1 h -1 ) were highest at Station 8 in July and 409 lowest at Station 11 in April (Fig 4; Table 4). The values of α β and P β max were positively correlated with each 410 other (r 2 = 0.99, n=6, Supplement 2).

412
The calculated production rate (

424
The daily benthic production rates (mg C m -2 d -1 ) were highest (> 495 mg C m -2 d -1 ) on 24 July 2019 in station 8 425 and lowest (< 1 mg C m -2 d -1 ) on 11 April 2019 for station 11 (Table 5). Daily production rates varied for 426 different assumptions with respect to vertical distribution of benthic chlorophyll-a in the sediment, with model 1 427 always supplying the highest values followed by model 2 (88%  10% compared to model 1) and model 3 (62%

428
 18% compared to model 1) (Table 5). The largest relative difference between models was found on 18

429
September for station 14, with the values for model 3 being 58% lower than those for model 1 (Table 5).

430
No significant relationship was found between daily production rates (mg C m -2 d -1 ) and corrected chlorophyll-a 431 concentrations (mg m -2 ) for any of the vertical distribution models, with the variation in biomass of benthic 432 algae explaining 23% (model 3) to 27% (model 2) of the variation in production rates (Supplement 3).   (Table 7). The values of the intercepts (65  13; 54  27) and slopes (118  49; 166  81) of these 453 relationships were comparable for these two sampling periods (Table 7).

454
For September 2018, no such relationship was found. For this reason, to determine the NDVI to CHLa 455 relationship that will be applied to satellite data only data collected in April and July 2019 will be used.

456
Relationships of corrected benthic chlorophyll-a concentrations as a function of NDVI and showed that the 457 values of NDVI as determined by means of hyperspectral sensors (based upon reflectance at 675 and 750 nm) were higher than the NDVIs as based upon red and near-infrared (NIR) spectral bands of the Landsat 7 ETM, 459 the Landsat 8 OLCI and the Sentinel 2 (Table 7).

460
For all three satellites under consideration, the relationships of the sum of corrected chlorophyll-a and 461 pheophytin-a concentrations as a function of NDVI explained more of the variance (56%-57%) than those of 462 relationships between corrected chlorophyll-a (44%) and uncorrected chlorophyll-a (47%) concentrations as a 463 function of NDVI (

487
Average concentrations were relatively high (more than 120 mg (CHLa+PHEOa) m -2 ) in late winter (27 488 February 2019) (Fig 5 & 6). From this time series, it appears that the highest peak on 21 April (2018) was 489 followed by a rapid decline in average concentrations of more than 100 mg (CHLa+PHEOa) m -2 in less than 50 490 days (Fig 6).

491
Concentrations of microphytobenthos pigments were, on average, relatively high at the lower parts (below NAP

512
Relationships between NDVI and pigment concentrations 513 No significant relationship between chlorophyll-a concentrations (mg m -2 ) and NDVI was found for the 514 campaign in September 2018, as was found for April and July 2019 (Table 7). This could have been caused by a 515 significant time lag (10-25 min) between the pigment sampling and NDVI measurements in September 2018.

516
[13] performed a controlled experiment in which artificial light illuminating a microphytobenthos layer was 517 switched on an off while NDVI was continuously monitored. A sharp increase in biomass after switching on the  (Table 7). The values of these slopes are relatively low compared to the range (from 229 to 1129 mg 530 m -2 ) and average values (523  236 mg m -2 ) as found in other studies (Table 1).

531
The highest correlation between NDVI and CHLa was obtained when CHLa was corrected (CHLa_c) and

546
In December 2018 and 2019, daily irradiances at noon fell below Ek (248 ± 164 µE PAR m -2 s -1 ; Table 6) so 547 light conditions were most probably too low to sustain net MPB growth (Fig 7). In February 2018, ice masses 548 covered the tidal flats which may have limited light penetration, although MPB biomass (this paper) and growth 549 has been observed under such circumstances [77]. In February 2019, when no ice was present, MPB biomass at 550 the higher tidal flats was higher than the year before (Fig 6). Biomass of MPB peaked in April, which is in line 551 by findings by [19]. From March to September, daily irradiances at noon were higher than Eopt (1915 ± 226 µE 552 PAR m -2 s -1 ; Table 6), implying that growth would be limited due to the light conditions if microphytobenthos 553 would then be at the sediment surface of emerged tidal flats. In January, February, October and November, 554 growth would be stimulated during relative sunny days (Fig 7).

570
[82], the daily maximum air temperature data of 2018 and 2019 could potentially explain the summer dip by 571 thermo-inhibition of MPB growth (Fig 7 & 8

581
Wind stress is generally higher during winter than during summer (Fig 7). Assuming a linear relationship  The consistent seasonal pattern as has been observed for the Dollard in this study underlines the feasibility of 598 using satellite images with a frequency of 5 to 10 images per year to monitor seasonal and, subsequently, year-

599
to-year variation in biomass of microphytobenthos on tidal flats of temperate coastal ecosystems. In particular, 600 if a long-term time series is being developed that enables decomposition of long-term trends and seasonal 601 dynamics.

602
Spatial patterns in NDVI

603
In the current study at Heringsplaat, highest concentrations of microphytobenthos were generally found at the 604 higher tidal flats (Fig 5-6  The EP-model allows for an estimate of a photo-inhibition parameter, but in the current study, photo-inhibition 618 did not occur (Fig 4). The maximum production rate is estimated to occur at a light intensity between 1493 and 619 2000 µE m -2 s -1 PAR ( 638 estimated in the current study ranged from 4.6 to 6.4 mm -1 (Table 4). This means that 90% of the light was lost 639 within the first 500 µm of the sediment.

640
The modelled distribution in this layer will thus determine the potential production rate in each sediment as it attenuation is estimated to be high, but production under these conditions can still be relatively high.
In the current study, however, the correlation between mud content and chlorophyll-a concentration was poor to 654 negative (Fig 3). Using model 3, under conditions with a high chlorophyll-a concentration and a low mud 655 content, a low mud content should result in a more homogenous distribution of chlorophyll-a in the sediment, 656 but with a high attenuation coefficient due to the high chlorophyll-a concentration. On these occasions, the daily 657 production is most likely to be under-estimated.

658
Temperature 659 The daily production rate for April 2019 is lower than can be expected based on the chlorophyll-a concentration 660 in the sediment (Table 5). A possible explanation might be found in the low temperatures at that date (Table 5).

661
In the current study, samples were incubated to estimate uptake rates of 14 C using in-situ water temperatures.

662
Production rates were then calculated using these water temperatures as well.

671
For the current study, the above would imply that in the April daily production rate is underestimated; while the 672 in-situ water temperature was low (7.6°C) actual sediment temperatures likely increased during exposure 673 becoming closer to the optimum temperature for production. In July, the opposite might have been the case, the 674 in-situ water temperature was 20°C and a decrease in production rate is expected during low tide due to sub-675 optimal sediment temperatures and the daily rates might have been over-estimated.

677
Matched optical chlorophyll and NDVI data from our field surveys showed that a relationship between both 678 indexes is feasible, although more campaigns would be desirable to increase the robustness and confidence. The

679
fact that the slope of the regression (Table 7) was lower than found by others (

692
Since NDVI can also be estimated from satellite data, our results imply that it should be possible to use satellite 693 images to estimate benthic primary production of relatively small and turbid estuaries. The Landsat sensors 694 monitor at a spatial resolution of 30m and have been providing data since the mid 1980's at bi-weekly 695 periodicity. Today, the European Sentinel-2 sensors as well as the ongoing Landsat missions 7 and 8, as well as 696 the forthcoming Landsat 9, provide a very valuable dataset that has yet to be analysed and interpreted. Based upon the three sampling campaigns, the average biomass of microphytobenthos was 131 ± 45 mg CHLa 699 m -2 (spectrophotometric; Table 3). Based upon the satellite images, annually averaged microphytobenthic 700 biomass in the Dollard was approximately 115 mg CHLa+PHEOa m -2 (Fig 5-6), and 86 mg CHLa m -2 under the 701 assumption that 75% of these pigments was CHLa (see Table 3

716
Based upon the three sampling campaigns, the average primary production of microphytobenthos ranged 717 between 73 mg C m -2 d -1 and 108 mg C m -2 d -1 , depending on the assumption on the vertical distribution of 718 chlorophyll-a in the sediment (Table 5). If assumed that these values represent averages for the growing season 719 and that the growing season runs from mid-January to mid-November (10 months), then the annual production 720 of microphytobenthos in the Dollard would have been between 22 g C m -2 y -1 and 33 g C m -2 y -1 . In 1976-1978,

721
the annual primary production of microphytobenthos of three stations in the Dollard was between 127 g C m -2 y -  With respect to seasonal scope for growth of microphytobenthos in temperate turbid environments as our study 735 area, insufficient light only restricts growth the benthic microalgae in mid-winter whilst growth may be 736 suboptimal due to photoinhibition in mid-summer (Fig 8). Under the assumption that MPB does not grow when 737 temperatures fall below 0C or rise above the maximum temperature tolerance, then growth is restricted in late 738 winter and in late summer. Wind stress was generally high from October to March and low in from April to Sept, whilst rainy days possibly disrupting the microphytobenthic layer occurred occasionally throughout the 740 year (Fig 7; Supplement 8).

741
Our findings suggest that optimal conditions for blooms of microphytobenthos in temperate turbid estuaries 742 occur in February or March, starting after sub-zero temperatures in (late) winter and depending on the 743 occurrence of storms and extreme rainfalls (Fig 7 & 8). Later in spring, both insolation and temperatures 744 become suboptimal for MPB growth, and growth might even become restricted from June to August. When 745 light and temperatures are declining after summer, an autumn bloom might still be possible in August and

746
September, but this cannot be as high as the spring bloom due to the ongoing nitrogen limitation and the onset 747 of the storm season (Fig 7 & 8; Supplement 8).

748
Seasonal dynamics in intertidal environments in north-western Europe are projected to change, e.g., due to a 749 decrease in the frequency of severe winters, an increase in the frequency of heat waves, more variations in 750 salinity due to extreme rainfall events in summer and, possibly, a decrease in ambient light conditions due to sea 751 level rise [93]. Subsequently, biomass and production of MPB is also likely to change, for example by an 752 advanced spring bloom and further growth restrictions during summer, resulting in a shift from unimodal [19] to 753 bimodal seasonality in MPB biomass and production (this paper). These changes will not only lead to changes 754 in spatiotemporal patterns of benthic primary production but also to changes in biodiversity of life under water 755 and ecosystem services including food supply.

757
During the drafting of the manuscript, Dr. Jacco C. Kromkamp

763
Evaline van Weerlee also performed sediment and chlorophyll-a analysis in the laboratory. Karel Bakker 764 measured DIC concentrations. Sonja van Leeuwen assisted with modelling of the primary productivity.

765
WaterProof B.V. and the Fieldwork Company provided transport to the tidal flats. We are most grateful to Peter

766
Herman and Elisabeth Addink for constructive comments on earlier versions of this manuscript.