Mapping the Intertidal Microphytobenthos Gross Primary Production, Part II: Merging Remote Sensing and Physical-Biological Coupled Modeling

Microphytobenthos (MPB) at the sediment surface of intertidal mudflats are known to show a high spatial and temporal variability in response to the biotic and abiotic conditions prevailing at the mud surface. It makes long-term and large-scale monitoring of MPB Gross Primary Production (GPP) difficult to set up. In this study, we developed the first 3D physical-biological coupled model (MARS-3D) that explicitly simulates GPP of intertidal MPB at the mudflat scale, and we compared the outputs with in situ and space remote sensing MPB GPP data. We discuss the sources of discrepancies between the modeling and the remote sensing approach in the light of future developments to be done. For instance, the remote sensing algorithm provides a very synoptic view of the mudflat GPP. It is well-suited to achieve diagnostic estimates of MPB GPP at the synoptic spatial and temporal scale. By contrast, the MARS-3D model provides a more dynamic representation of the MPB activity and prognostic estimates of MPB GPP over the mudflat. It is very relevant to resolve the seasonal and inter-annual dynamics of MPB. Getting comparable GPP estimates derived from the remote sensing algorithm and 3D physical-biological coupled model will further require a better convergence in terms of equations structure, biological constants parameterization, and source data used (i.e., vegetation index vs. chlorophyll a). Setting a common parameterization in both the numerical model and remote sensing algorithm might be challenging in a perspective of mapping MPB PP over large mudflats from a synoptic to inter-annual time scale, but it could open the door to a new way of quantifying MPB GPP over large intertidal mudflats.

Microphytobenthos (MPB) at the sediment surface of intertidal mudflats are known to show a high spatial and temporal variability in response to the biotic and abiotic conditions prevailing at the mud surface. It makes long-term and large-scale monitoring of MPB Gross Primary Production (GPP) difficult to set up. In this study, we developed the first 3D physical-biological coupled model (MARS-3D) that explicitly simulates GPP of intertidal MPB at the mudflat scale, and we compared the outputs with in situ and space remote sensing MPB GPP data. We discuss the sources of discrepancies between the modeling and the remote sensing approach in the light of future developments to be done. For instance, the remote sensing algorithm provides a very synoptic view of the mudflat GPP. It is well-suited to achieve diagnostic estimates of MPB GPP at the synoptic spatial and temporal scale. By contrast, the MARS-3D model provides a more dynamic representation of the MPB activity and prognostic estimates of MPB GPP over the mudflat. It is very relevant to resolve the seasonal and inter-annual dynamics of MPB. Getting comparable GPP estimates derived from the remote sensing algorithm and 3D physical-biological coupled model will further require a better convergence in terms of equations structure, biological constants parameterization, and source data used (i.e., vegetation index vs. chlorophyll a). Setting a common parameterization in both the numerical model and remote sensing algorithm might be challenging in a perspective of mapping MPB PP over large mudflats from a synoptic to inter-annual time scale, but it could open the door to a new way of quantifying MPB GPP over large intertidal mudflats.

INTRODUCTION
Benthic microalgae or microphytobenthos (MPB) inhabiting the sediment surface sustain a high biological production in intertidal mudflats (MacIntyre et al., 1996;Underwood and Kromkamp, 1999). As the main primary producers on intertidal mudflats, MPB are of key importance for higher trophic levels from benthic fauna to birds (Herman et al., 2000;Kang et al., 2006;Jardine et al., 2015) and for pelagic organisms that benefit MPB locally resuspended by tides and waves (Perissinotto et al., 2003;Krumme et al., 2008) but also exported to adjacent ecosystems (Saint-Béat et al., 2013). With a global annual Primary Production (PP) estimated to ∼500 million tons of Carbon (C; Cahoon, 1999), MPB also participate in the Blue Carbon (Otani and Endo, 2019). Guarini et al. (2008) suggested that MPB PP could represent a significant amount of carbon not considered in the global carbon cycle. However, their contribution to the global carbon budget remains unknown.
The spatial and temporal distribution of MPB over mudflats is highly variable, as it is driven by highly variable physical [light, mud surface temperature (MST), tides, waves, and current] and biological (grazing, biostabilization, and bioturbation) conditions (e.g., Admiraal, 1984;Blanchard et al., 1996;MacIntyre et al., 1996;Underwood, 2001;Morris and Kromkamp, 2003;Sahan et al., 2007;Salleh and McMinn, 2011;Kwon et al., 2014;Orvain et al., 2014a,b;Savelli et al., 2019). Such a variability impedes any accurate and robust assessment of the role played by MPB at the scale of the whole mudflat ecosystem and of its contribution to the carbon cycle. MPB PP and biomass measurements are usually limited to single-point sampling (e.g., Vieira et al., 2013;Orvain et al., 2014a;Cartaxana et al., 2015;Pniewski et al., 2015). This approach succeeds in capturing the MPB temporal dynamics but is rapidly limited when dealing with spatial and temporal variations of MPB PP and biomass at the scale of an entire mudflat. Only a few studies resolved the MPB spatial variability at the mudflat scale, as time and important logistical resources are required to meet this goal (e.g., Guarini et al., 1998;Ubertini et al., 2012). Remote sensing and physical-biological coupled modeling are then relevant and non-invasive approaches to infer on MPB dynamics (e.g., Guarini et al., 2000;Combe et al., 2005;van der Wal et al., 2010;Brito et al., 2013). Jobson et al. (1980) initiated the use of remote sensing to assess the MPB biomass from a tower-mounted sensor designed to scan a mudflat of South Carolina (USA). Since then, airborne and space remote sensing methods were increasingly developed and more widely used in MPB studies (e.g., Méléder et al., 2003a;Brito et al., 2013;Benyoucef et al., 2014;Daggers et al., 2018). Remote sensing data can cover large spatial scales (∼ from one meter to several kilometers), and multispectral broadband sensors promise high quality data to map MPB biomass and PP over entire mudflats. Daggers et al. (2018) first coupled in situ measurements, satellite remote sensing data, and observed tidal heights to map synoptic MPB PP in spring at the scale of the Oosterschelde and Westerschelde estuaries (The Netherlands). Méléder et al. (2020) coupled in situ measurements, satellite remote sensing data, and data simulated (light, MST, and tidal height) by a 3D physical-biological coupled model to map synoptic MPB PP at three seasons over the large intertidal Brouage mudflat on the French Atlantic Coast.
Recently, remotely sensed estimates of vegetation index and of in water MPB chlorophyll a (Chl a) concentration were compared to model outputs in order to assess the model ability to simulate realistic MPB biomass levels over the Brouage mudflat (Savelli et al., 2018(Savelli et al., , 2019. Such a comparison does not exist for MPB PP. The recent development for the Brouage mudflat of a regional MPB Gross Primary Production (GPP) algorithm (Méléder et al., 2020) and of a 1D MPB GPP physical-biological coupled model (Savelli et al., 2018) allows for the comparison of MPB GPP estimates derived from space remote sensing and a prognostic modeling approach. The objective of this study is to infer the capacity of remote sensing and prognostic modeling to converge toward realistic MPB GPP estimates over the large Brouage mudflat. In this paper, we describe first the physicalbiological coupled model developed in 3D in order to simulate the spatial and temporal dynamics of intertidal MPB. Second, we compare the model outputs with remotely sensed MPB GPP estimates coincident in space and time. Finally, we discuss the sources of discrepancies between the two approaches in the light of future developments to be done.

Study Site
The Pertuis Charentais Sea is a shallow semi-enclosed sea where develops one of the biggest shellfish farming activity in Europe (Goulletquer, 1998). It receives riverine inputs originating from the agricultural watersheds of the Sèvre, Charente, and Seudre rivers (Figure 1). Located in the southern part of the Pertuis Charentais, the Brouage intertidal mudflat is a 42-km 2 intertidal mudflat composed of fine cohesive sediments (Bocher et al., 2007) distributed over a gentle slope (∼ 1/1,000; Le Hir et al., 2000). As in many other mudflats along the northern European Atlantic coast, a dense MPB biofilm develops at the sediment surface at low tide with Chl a concentrations reaching up to 25 mg Chl a m −2 (Guarini, 1998). The Brouage mudflat is responsible for a large part of the high PP reported in the Marennes-Oléron Bay (Struski and Bacher, 2006).

Observations
Field campaigns were conducted during spring and daytime low tides on 5-6 May and 2-3 July 2015. A total of 9 GPP estimates were derived from carbon fluxes at the air-sediment interface measured with benthic chambers. The CO 2 concentration was measured in the chambers continuously over a period of 20 to 30 min using an infrared gas analyzer (IRGA). At the same time, incident photosynthetically active radiation (400 to 700 nm; PAR, µmol photons m −2 s −1 ) and temperature ( • C) were measured (LS-C planar sensor plugged to a ULM-500 data-logger, Walz, Effeltrich, Germany) at the sediment surface near the chambers at a 30-s frequency. The MPB biomass was estimated continuously (every 5-10 min) by sediment sampling of the upper 250 µm layer by the "crème brulée" method (Laviale et al., 2015). 43 biomass samples were taken from the sediment surface. The Chl a content of sediment was determined by reversed phase HPLC (Hitachi High Technologies Co., Japan) calibrated for Chl a. The Chl a concentration was normalized to the sampled surface (1.5 cm 2 ) to be expressed in mg Chl a m −2 . The sampling protocol is fully detailed in Méléder et al. (2020).

The MARS-3D Model
The MARS-3D model (3D hydrodynamical Model for Applications at Regional Scale) is a regional ocean model that simulates the ocean physics (Lazure and Dumas, 2008). In this study, we used the regional configuration set up for the Pertuis Charentais area, including the Brouage mudflat. The model was discretized into 100 m by 100 m horizontal grid cells and 20 sigma-levels over depth. The model was run for the same domain as in Polsenaere et al. (2017) (Figure 1). The MARS-3D model is fully detailed in Lazure and Dumas (2008). Atmospheric forcings (10 m winds, air temperature, atmospheric pressure at sea level, nebulosity fraction, relative humidity, and downwelling solar fluxes) were provided by the Meteo France AROME model (https://donneespubliques.meteofrance.fr). At the open boundaries of the numerical grid, the MARS-3D model was forced by tidal amplitudes and phases of 115 harmonic constituents from the cstFRANCE tidal model developed by the French marine service for hydrography and oceanography (SHOM; Simon and Gonella, 2007). Initial and boundary conditions of seawater temperature, salinity, current velocity, and sea surface height were derived from the MANGAE 2500 Ifremer model (Lazure et al., 2009).

The Mud Surface Temperature Model
The mud surface temperature model used in Savelli et al. (2018) was coupled to MARS-3D. Thermodynamic equations detailed in Savelli et al. (2018) simulated heat fluxes within a 1 cm deep sediment layer. No horizontal fluxes were considered. In their study, Savelli et al. (2018) successfully compared the simulated MST with 1 min MST data measured in situ on the Brouage mudflat. The differential equation of heat energy balance was solved by the MARS-3D numerical scheme. The MST model is fully detailed in Savelli et al. (2018).

The MPB Model
The MPB model used in Savelli et al. (2018) was also coupled to MARS-3D. The MPB model simulated the MPB biomass in both the surface biofilm (S, mg Chl a m −2 ) and sediment first centimeter (F, mg Chl a m −2 ), and the grazer biomass (Peringia ulvae, Z, and mg C m −2 ) at the sediment surface. The MPB model accounted for vertical MPB migrations driven by diurnal and tidal cycles through exchanges of MPB biomass between S and F (Guarini et al., 2000). The P. ulvae growth was sustained by grazing on the MPB biofilm. It was controlled by the MST and MPB biomass in the biofilm. The biomass-specific photosynthetic rate P b [mg C (mg Chl a) −1 h −1 ] was regulated by MST and PAR according to the models of Blanchard et al. (1996) and Platt and Jassby (1976), respectively. In the present study, the MPB biomass in the biofilm referred to the variable S * introduced by Savelli et al. (2018) that represented the S compartment that incorporated the S instantaneous production of biomass [mg Chl a m −2 ], which is directly transferred to F. The MPB model and differential equations are fully detailed in Savelli et al. (2018).
The physical-biological coupled model was initialized with a spin-up starting from 12 September 2014 00:00:00 UTC to 1 January 2015 00:00:00 UTC. The variables F, S, and Z were initially set to 100 mg Chl a m −2 , 0 mg Chl a m −2 and 1,000 mg C m −2 , respectively. The physical-biological coupled model was then run from 1 January 2015 00:00:00 UTC to 1 January 2016 00:00:00 UTC.

Remotely Sensed MPB GPP
The MPB GPP algorithm developed by Méléder et al. (2020;GPP-algo) coupled Normalized Difference Vegetation Index data (NDVI; Tucker, 1979) derived from the SPOT and Pléiades satellite sensors with MARS-3D derived forcings. In the GPPalgo, the remotely sensed GPP was obtained by constraining with tidal heights, light levels and mud surface temperature simulated by MARS-3D the horizontal distribution of MPB biomass estimated from the NDVI data and modulated by MPB vertical migration. The photosynthetic rate of MPB in the GPP-algo [P b , mg C (NDVI) −1 m −2 h −1 ] was estimated by the temperature and light-related production models of Blanchard et al. (1996) and Eilers and Peeters (1988) parameterized with photophysiological parameters fitted on the laboratory measured light curves (α the initial slope of the curve, E opt the optimum irradiance for photosynthesis, and P b max the photosynthetic capacity). The GPPalgo is fully detailed in Méléder et al. (2020).

Comparison of Remotely Sensed and Simulated MPB GPP
We compared the MPB GPP remotely sensed and simulated by MARS-3D on satellite acquisition matching days in May and July 2015. We ran two model setups. In the first run, the maximum photosynthetic capacity [P b MAX , mg C (mg Chl a) −1 h −1 ] was seasonally estimated by Blanchard et al. (1997; MARS-3D season run). In a second run, P b MAX was set from in situ GPP measurements in May and July 2015 (MARS-3D synoptic run). P b MAX was estimated from the measured biomassspecific production rate P b [mg C (mg Chl a) −1 h −1 ], light, and MST from which we retrieved values of P b MAX with the models of Blanchard et al. (1996) and Platt and Jassby (1976) parameterized with photophysiological parameters (β the shape parameter of the production-temperature relationship, T opt the temperature optimum for MPB photosynthesis, T max the maximum temperature for MPB photosynthesis and E k the light saturation parameter) from Savelli et al. (2018). We iterated values of P b MAX by a dichotomous analysis based on the intermediate value theorem (Bolzano, 1817). In MARS-3D synoptic , the MPB GPP was simulated in days matching the in situ measurements with MARS-3D parameterized with the mean value of P b MAX iterated during in situ measurements in May and July 2015.
First, we took advantage of the MARS-3D season model to investigate the seasonal dynamics of MPB biomass and GPP (Figure 2). Second, we explored the spatial distribution of MPB biomass and GPP simulated by MARS-3D season (Figure 2). Then, we compared, in a single-point approach, P b MAX and GPP extracted from the MARS-3D season and MARS-3D synoptic grid cell corresponding to the sampling site with in situ measurements and GPP-algo (Figure 2). Finally, we compared fields of MPB GPP obtained with MARS-3D synoptic and GPP-algo (Figure 2).

Seasonal MPB Dynamics Simulated at the Sampling Site
In the MARS-3D season run, the MPB biomass simulated in the sediment 1 st cm reached one seasonal maximum on 25 February and 30 December with ∼190 mg Chl a m −2 in 2015 ( Figure 4A). The seasonal minimum of MPB biomass simulated in the sediment occurred on 22 August with 34.2 mg Chl a m −2 ( Figure 4A). The MPB biomass in the biofilm simulated in the MARS-3D season run was 20.6 ± 11.25 mg Chl a m −2 at the study site and varied from 0 to 44.3 mg Chl a m −2 ( Figure 4B). The mean hourly GPP during daytime emersion simulated in the MARS-3D season run was 71.3 ± 65.4 mg C m −2 h −1 ( Figure 4C and Table 1). Similarly to MPB biomass simulated in the biofilm, it was highly variable ranging from 0 to 278.8 mg C m −2 h −1 and often reached GPP levels similar to those measured in May and July 2015 ( Figure 4C). The mean daily GPP simulated in the MARS-3D season run at the study site was 359.9 ± 229.5 mg C m −2 d −1 (Table 1). Such simulated GPP rates resulted in an annual GPP at the study site of 131 g C m −2 yr −1 ( Table 1).

MPB GPP Simulated Over the Whole Mudflat
The MPB biomass in the biofilm simulated in MARS-3D season on satellite acquisition matching days in May 2015 was ∼40 and 35 mg Chl a m −2 on the upper and lower shore, respectively ( Figure 5A). The MPB biomass simulated in the biofilm on satellite acquisition matching days in July 2015 was relatively homogeneous over the mudflat and was lower than in May 2015 with ∼25 and 20 mg Chl a m −2 ( Figure 5B). In May 2015, the hourly GPP simulated in the MARS-3D season run on satellite acquisition matching days was relatively homogeneous over the mudflat with values of 130 mg C m −2 h −1 (Figure 6A). In July 2015, the hourly GPP simulated in the MARS-3D season run on satellite acquisition matching days was higher on the southern part of the mudflat (∼ 40 mg C m −2 h −1 ) than on the northern part (∼ 20 mg C m −2 h −1 ; Figure 6B). In May 2015, the daily integrated GPP simulated in MARS-3D season was higher on the upper shore (∼900 mg C m −2 d −1 ) than on the lower shore (∼700 mg C m −2 d −1 ; Figure 6C). Integrated over the mudflat, GPP simulated in MARS-3D season during satellite acquisition matching day in May 2015 was 19.5 t C (Figure 6C). In July 2015, similarly to the hourly GPP, the daily GPP simulated in the MARS-3D season run on satellite acquisition matching days was higher on the southern part of the mudflat (∼ 250 mg C m −2 d −1 ) than on the northern part (∼ 100 mg C m −2 d −1 ; (Figure 6C). It represented 3.8 t C, once integrated over the mudflat (Figure 6D).

MPB GPP Single-Point Comparison: Simulated vs. Remotely Sensed and in situ Data
MAX retrieved from iterations on in situ measurements was in average 0.26 ± 0.11 and 0.67 ± 0.30 mg C (mg Chl a) −1 h −1 in May and July 2015, respectively ( Figure 7A and Table 2). In the MARS-3D season run, P b MAX was ∼9.4 and 6.4 mg C (mg Chl a) −1 h −1 in days matching the in situ measurements in May and July 2015, respectively ( Figure 7A and Table 2). It was hence 36-and 9-fold higher than P b MAX retrieved from in situ measurements in May and July 2015, respectively ( Figure 7A and Table 2). In the MARS-3D season run, the MPB biomass simulated in the biofilm during the field campaign was on average 21.7 ± 12.5 and 21.1 ± 9.5 mg Chl a m −2 in May and July 2015, respectively ( Figure 7B). It was not significantly different than the measured MPB biomass in the biofilm in May and July 2015 (Mann-Whitney test: p = 0.31 and 0.95, respectively; Figure 7B). The simulated GPP in the MARS-3D season run in days matching the in situ measurements (164.8 ± 66.7 mg C m −2 h −1 ) was on average 29-to 40-fold higher than the measured (5.69 ± 3.22 mg C m −2 h −1 ) and remotely sensed (4.13 ± 2.22 mg C m −2 h −1 ) GPP in May 2015, respectively (Mann-Whitney test: p < 0.01; Figure 7C and Table 2). In July 2015, GPP simulated in the MARS-3D season run (41.3 ± 43.6 mg C m −2 h −1 ) was significantly different than in situ (6.3 ± 0.3 mg C m −2 h −1 ) and remotely sensed GPP (2.2 ± 1.4 mg C m −2 h −1 ; Mann-Whitney test: p < 0.01; Figure 7C and Table 2).
In the MARS-3D synoptic run, i.e., with a mean P b MAX based on in situ measurements in May and July 2015 (0.26 ± 0.11 and 0.67 ± 0.30 mg C (mg Chl a) −1 h −1 , respectively; Figure 8A), the MPB biomass simulated in the biofilm was consistent with the estimates measured in situ (Mann-Whitney test: p = 0.6 and 0.62, respectively; Figure 8B). The mean MPB biomass simulated in the MARS-3D synoptic run on in situ sampling days in May 2015 was 17.8 ± 9.9 mg Chl a m −2 (Figure 8B). In July 2015, it was 17.3 ± 9.9 mg Chl a m −2 (Figure 8B). In the MARS-3D synoptic run, the simulated GPP compared to GPP measured in situ and derived from GPP-algo on in situ sampling days. With values of 5.1 ± 2.13 mg C m −2 h −1 in May 2015, GPP simulated in the MARS-3D synoptic run was not significantly different than GPP measured in situ and derived from GPP-algo (Mann-Whitney test: p = 0.9 and 0.5; Figure 8C and Table 2). In July 2015, GPP simulated in the MARS-3D synoptic run (5.25 ± 4.78 mg C m −2 h −1 ) was not significantly different than GPP measured in situ and derived from GPP-algo (Mann-Whitney test: p = 0.6 and 0.08, respectively in July 2015; Figure 8C and Table 2).

Variability
The parametrization of P b MAX in the MARS-3D synoptic run resulted in much lower simulated GPP values over the whole mudflat on satellite acquisition matching days in May and July 2015 than in the MARS-3D season run (Figure 9). In May 2015, the daily integrated GPP simulated in MARS-3D synoptic was homogeneous over the mudflat reaching value of ∼22 mg C m −2 d −1 (Figure 9A). Integrated over the whole mudflat, GPP was 0.54 t C in the MARS-3D synoptic run. In July 2015, the daily integrated GPP simulated in the MARS-3D synoptic run was higher on the upper shore (∼24 mg C m −2 d −1 ) than on the lower shore (∼15 mg C m −2 d −1 ; Figure 9B). It resulted in a simulated spatially integrated GPP of 0.40 t C in the MARS-3D synoptic run ( Figure 9A).
Compared to the remotely sensed GPP data, the use of P b MAX based on synoptic field data in the MARS-3D synoptic run resulted into slightly lower daily integrated GPP simulated over the mudflat in May 2015 (−21.2 ± 58.2 mg C m −2 d −1 ; Figure 9E). Integrated over the mudflat, the simulated GPP decreased by 0.88 t C compared to the remotely sensed GPP estimate ( Figure 9E). Aside from the extreme upper shore of the southern part of

Variables Units Values
Hourly GPP mg C m −2 h −1 71.3 ± 65.4 Daily GPP mg C m −2 d −1 359.9 ± 229.5 Annual GPP g C m −2 yr −1 131 the mudflat where MARS-3D did not simulate MPB, the GPP differences were particularly high on patches in the central and northern part of the mudflat with differences up to 40 mg C m −2 d −1 in May 2015 ( Figure 9E). The GPP differences between the MARS-3D synoptic run and the GPP-algo were lower in July 2015 than in May 2015. The mean GPP difference was −3.82 ± 15.7 mg C m −2 d −1 in July 2015 ( Figure 9F). The MARS-3D synoptic /GPP-algo difference of mudflat-integrated GPP in July 2015 was −0.40 t C ( Figure 9F). The high temporal variability of P b MAX impeded therefore a convergence of its estimation between the modeling and space remote sensing approach. Consequently, it was a very sensitive parameter in the model as it mediated strong differences in GPP estimates between the MARS-3D synoptic and MARS-3D season runs.

Simulated and Remotely Sensed MPB GPP Estimates
The GPP simulated in the MARS-3D season run model is much higher than GPP measured in situ. Over the mudflat, the GPP simulated in the MARS-3D season run is also higher than the GPP derived from the GPP-algo developed by Méléder et al. (2020) for the whole mudflat. However, both the model (MARS-3D season ) and the remote sensing algorithm provide hourly and daily GPP rates in the range of GPP values reported in the literature (Cahoon, 1999;Underwood and Kromkamp, 1999;Daggers et al., 2018).
Given that the MARS-3D model is constrained by simulated water height and meteorological parameters, it is sensitive to likely inaccuracies in the forcings that might impede the model to resolve the high temporal variability of the physical environment. Nevertheless, simulated PAR and MST data lie within the range of in situ measurements and the impact of such inaccuracies on GPP estimates may be limited. As the MPB biomass simulated in the biofilm in the MARS-3D season run is also consistent with in situ measurements, the MARS-3D season run-observations GPP discrepancies can be attributed to differences in the MPB maximum photosynthetic capacity (P b MAX ) used in the   MARS-3D season run and estimated in the field in May and July 2015 at the study site. When set up with P b MAX comparable to the measured values, the GPP simulated in the MARS-3D synoptic run better compared to in situ GPP measurements. This suggests that, at the mudflat scale, MPB GPP estimates derived from remote sensing and the model are sensitive to the MPB photophysiological parameters (P b MAX but also the other temperature and light-related photosynthesis parameters such as the temperature optimum and maximum for MPB photosynthesis, the shape parameter of the production-temperature relationship and the light saturation parameter) and their spatio-temporal variability.

From Synoptic to Seasonal GPP Estimates
The GPP-algo developed by Méléder et al. (2020) is parameterized with synoptic measurements of the photosynthetic activity and the related photophysiological parameters (P b MAX and also the optimal irradiance for photosynthesis and the initial slope of the production-irradiance relationship) of MPB cells collected during the field campaigns. It is therefore well-suited to depict the high temporal variability of the MPB photosynthetic response to the physical environment. However, GPP estimates from space remote sensing are restricted to the satellite data availability, which depends on the satellite revisit time, the cloud cover and the time window of acquisition during the day (Daggers et al., 2018;Méléder et al., 2020). Despite this limitation, remote sensing GPP algorithms are relevant to estimate MPB GPP at the synoptic time scale (Daggers et al., 2018;Méléder et al., 2020).
When parameterized with synoptic in situ estimate of P b MAX , the MARS-3D model simulates GPP values that also compare to the in situ estimates. However, such a parametrization does only apply to a specific location at a specific time. Despite GPP simulated in MARS-3D season on in situ sampling days depart from measured GPP, daily and annual GPP simulated at the study site in 2015 (359.9 ± 229.5 mg C m −2 d −1 and 131 g C m −2 yr −1 , respectively) are consistent with the literature (Cahoon, 1999;Underwood and Kromkamp, 1999;Savelli et al., 2018Savelli et al., , 2019. Annual MPB GPP estimates can be obtained from extrapolation of daily GPP (3.65 to 93.99 g C m −2 yr −1 ; Méléder et al., 2020) derived from GPP-algo. However, the high MPB GPP variability at the hourly scale makes such extrapolations to be considered with caution. In contrast, the relatively consistent GPP values simulated at high frequency (12 s time step) in MARS-3D season over a year are likely to be used with more confidence for estimating GPP at the seasonal scale. The MPB model used in this study is adapted from the 1D model developed and validated in Savelli et al. (2018), which reasonably simulates the MPB dynamics in the Brouage mudflat for the year 2008. Similarly to Savelli et al. (2018), the seasonal cycle of MPB biomass simulated at the study site in 2015 is characterized by a spring bloom, a summer depression and a fall bloom. The fair agreement between the MPB biomass in the biofilm simulated in MARS-3D with the time-coincident observations suggests that overall the model simulates with some confidence the MPB dynamics at the seasonal scale in 2015.

From Single-Point to Mudflat GPP Estimates
MPB GPP estimates derived from remote sensing algorithms and physical-biological coupled models depend on the photophysiological parameters values and as such, on their sampling location on the mudflat. On a sandflat of the Bay of Paranaguá (Brazil), Fonseca et al. (2008) measured with benthic chambers higher PP rates in the upper and middle shores (1.9-2.1 g C m −2 d −1 and 1.3-2.2 g C m −2 d −1 , respectively) than in the lower shore (0.24-0.27 g C m −2 d −1 ). Cook et al. (2004) measured CO 2 fluxes at the air-sediment interface at two tidal levels of a mudflat located in Tasmania. The uptake of inorganic carbon (total CO 2 ) at the benthic interface was higher on the upper shore (up to 15,000 µmol m −2 h −1 ) than on the lower shore (up to 6,000 µmol m −2 h −1 ), suggesting a higher benthic GPP on the upper shore than on the lower shore. On the Brouage mudflat, the relatively low GPP over the whole area depicted by the GPP-algo (Méléder et al., 2020) may be the result of a parametrization of the algorithm based on photophysiological parameters estimated on potentially low-productive MPB cells collected on the lower shore. Conversely, the photophysiological parameters used in the MARS-3D parameterization were derived from MPB cells collected on the middle shore of the Brouage mudflat (Figure 1; Blanchard et al., 1997). Consequently, when applied to the entire mudflat, such a parametrization may result in a GPP overestimation, especially on the lower shore as suggested by the MARS-3D season -in situ measurements mismatch reported in the lower shore. Such a model-data mismatch is reduced in the MARS-3D synoptic run when the MARS-3D model is parameterized using MPB photophysiological parameters estimated on mud samples gathered on the lower Brouage shore.
Remote sensing algorithms and physical-biological coupled model GPP estimates rely on the photosynthetically active MPB biomass at the mud surface. While the remote sensing algorithm uses NDVI data, the MARS-3D model uses Chl a concentration simulated in the biofilm to infer on the horizontal distribution of MPB biomass at the mud surface. NDVI data provide a synoptic view of the MPB activity for a given time. Combining these NDVI snapshots over time (i.e., a diurnal cycle) requires us to account for the MPB vertical migration scheme during daytime low tides. Daggers et al. (2018) introduced the MPB vertical migrations in their remote sensing algorithm by modulating the PP rate during the first hour of the daytime emersion period. In their algorithm (GPP-algo), Méléder et al. (2020) assumed that the MPB biomass detected by satellite corresponds to the fully-established biofilm during the daytime low tide (total photosynthetically active biomass). Méléder et al. (2020) considered therefore a progressive establishment of the total photosynthetically active biomass at the sediment surface. In the MARS-3D model, the MPB biomass simulated in the biofilm follows the MPB vertical scheme described by Guarini et al. (2000). MPB cells migrate upward from the lower 1st cm sediment to the sediment surface during daytime low tides. At nightfall or at the time the flood begins, MPB cells migrate back downward. As the MPB biomass simulated in the biofilm in the MARS-3D model compares to the time-coincident field measurements, the model can resolve with some confidence the temporal variability of the MPB biomass in the biofilm. However, no gridded data of benthic Chl a are available to assess the ability of the MARS-3D model to resolve the spatial variability of the MPB biomass in the biofilm. In the MARS-3D model (MARS-3D season ), the MPB biomass simulated in the biofilm is slightly higher on the upper shore than on the lower shore on satellite acquisition matching days. However, the time-coincident NDVI data suggest a higher MPB biomass on the lower and middle shores than on the upper shore (Méléder et al., 2020). Comparing the MARS-3D model (Chl a) and GPP-algo (NDVI) is difficult as the NDVI-Chl a relationship is not linear, especially at high values of Chl a (Méléder et al., 2003a,b;Serôdio et al., 2009). The remote sensing of the MPB biomass in Chl a units from hyperspectral imaging might overcome the MPB biomass units mismatch (Kazemipour et al., 2012;Launeau et al., 2018).
The horizontal resolution of satellite data and 3D regional models is also a critical issue when estimating GPP of patchydistributed MPB. Méléder et al. (2020) report high differences in the NDVI signal between in situ and satellite observations and between different satellite sensors due to the dilution of the NDVI signal with the increasing pixel size and the patchiness distribution of the biofilm (Saburova et al., 1995;Spilmont et al., 2011). As the horizontal resolution of the MARS-3D model (100 m) and the remote sensing algorithm developed by Méléder et al. (2020) (from 2 to 6 m) differs, confronting quantitatively remotely sensed and simulated GPP per unit of surface must be considered with caution. The high horizontal resolution of remote sensing data is appropriate to monitor the MPB patchiness, which is not the case of the MARS-3D model. For this reason, remotely sensed GPP estimates are more suitable for comparison with synoptic in situ measurements.

CONCLUSIONS
This study is a first attempt to simulate the 3D MPB dynamics at the scale of an entire intertidal mudflat. Combined with a novel space remote sensing approach to assess MPB GPP, it allows for a first comparison of MPB GPP estimates derived from a remote sensing algorithm (GPP-algo) and a regional 3D physicalbiological coupled model. The remote sensing algorithm provides a very synoptic view of the mudflat GPP. It is well-suited to achieve diagnostic estimates of MPB GPP at the synoptic spatial and temporal scale. By contrast, the 3D physical-biological model provides a more dynamic representation of the MPB activity as well as prognostic estimates of MPB GPP over the mudflat. It is very relevant to resolve the seasonal and inter-annual dynamics of MPB. Furthermore, the coupling of the intertidal and pelagic domains in the regional 3D model could be envisaged in the future to assess the fate in the coastal ocean of fresh organic carbon resulting from MPB GPP. However, a refinement of its horizontal numerical mesh is required to resolve the MPB patchiness and to allow a better comparison with high resolution remote sensing data in the future. With respect to remote sensing, GPP algorithms are still limited by the too low spectral resolution of the multispectral (3-10 bands) satellite sensors and by the data availability. Hyperspectral remote sensing is able to capture photosynthetic capabilities and GPP, as recently proposed for terrestrial vegetation (DuBois et al., 2018;Lees et al., 2018). This approach starts to be develop successfully on MPB by Méléder et al. (2018). Furthermore, airborne hyperspectral data (hundreds of bands) could complement space satellite remote sensing data in an era of remote sensing drone aircraft democratization (Launeau et al., 2018). The start-up of the Deutsches Zentrum für Luft-und Raumfahrt Earth Sensing Imaging Spectrometer (DESIS) on-board of the International Space Station could also enable the development of Earth Observation algorithms based on hyperspectral images from space. Confronting GPP derived from remote sensing algorithms and 3D physicalbiological models will require a better convergence in terms of equations structure, biological constants parameterization, and source data used (i.e., NDVI vs. Chl a). Such a convergence would provide very complementary tools for diagnostic and prognostic analyses of the MPB GPP evolution at mudflat scales. While space remote sensing algorithms may provide a more realistic view of the MPB dynamics at the mudflat scale, 3D coupled physical-biological models can fill the gap left by space remote sensing strongly impacted at these latitudes by cloud cover, hence allowing for an annual budget of MPB GPP. Consequently, remote sensing algorithms and 3D coupled physical-biological models can be combined to monitor in an operational way MPB GPP from the synoptic to the annual scale and to achieve annual MPB GPP budget for large intertidal mudflats. Such a convergence was acclaimed in Babin et al. (2015) for phytoplankton in remote environments, whereas, for mudflats, remote sensing is pivotal for PP monitoring. Such an achievement will however require spatial and temporal surveys of the MPB photophysiological parameters across tidal heights in order to better assess the MPB photosynthetic response in time and space and better parameterize the remote sensing algorithms and models. Assessing the photosynthetic response of MPB to its highly variable environment is a challenge for the coming years in a perspective of quantifying MPB PP over large productive mudflats from a synoptic to inter-annual time scale.

DATA AVAILABILITY STATEMENT
Several chapters of RS.'s PhD thesis will use the model presented in this study. As a consequence, the model data presented in this study were archived in a ZENODO repository (https://zenodo. org/record/4022383), which will be available after an embargo period corresponding to the completion of this PhD thesis (December 2020).

AUTHOR CONTRIBUTIONS
RS, VM, VL, PC, PP, CD, and JL set the conceptual framework of this study. RS, PC, VL, and PP coupled the MPB model to MARS-3D. VM, RS, AB, PP, and JL carried out the experiments. RS, AB, and JL analyzed the data. RS wrote the manuscript. RS, VM, VL, PC, PP, CD, JL, and AB reviewed the manuscript. All authors contributed to the final version of the manuscript.

FUNDING
This work was supported by (i) the DYCOFEL project, funded through the 2015 Fondation de France call Quels littoraux pour demain?; (ii) the MIMOSA project, funded through the 2018 CNRS EC2CO-LEFE call; (iii) the HYPEDDY project, funded through the 2018-2020 Tosca-CNES call; (iv) the BIO-Tide project, funded through the 2015-2016 BiodivERsA COFUND call for research proposals, with the national funders BelSPO, FWO, ANR, and SNSF; (v) the public funds received in the framework of GEOSUD, a project (ANR-10-EQPX-20) of the program Investissements d'Avenir, managed by the French National Research Agency; and (vi) the projects Littoral 1 and ECONAT funded by the Contrat de Plan Etat-Région (CPER) and the CNRS and the European Regional Development Fund. Pléiades and SPOT images were acquired by CNES's ISIS program, facilitating scientific access to imagery. Pléiades CNES 2015, 2018, Distribution Airbus DS, all rights reserved. Commercial uses forbidden. This research was part of fulfillment of the requirements for a Ph.D. degree (RS) at the Université de La Rochelle, France. RS was supported by a Ph.D. fellowship from the French Ministry of Higher Education, Research and Innovation.