CDOM Spatiotemporal Variability in the Mediterranean Sea: A Modelling Study

: This study investigates the spatial and temporal variability of chromophoric-dissolved organic matter (CDOM) in the Mediterranean Sea. The analysis is carried out using a state-of-the-art 3D biogeochemical model. The model describes the plankton dynamics, the cycles of the most important limiting nutrients, and the particulate and dissolved pools of carbon. The source of CDOM is directly correlated to the dynamics of dissolved organic carbon (DOC) by a ﬁxed production quota. Then CDOM degrades by photobleaching and remineralization. The main innovation of the system is the inclusion of a bio-optical radiative transfer model that computes surface upwelling irradiance, and therefore simulates remotely sensed reﬂectance (Rrs). Simulation results of three model conﬁgurations are evaluated using satellite Rrs, particularly at 412 nm, 443 nm, and 490 nm. All simulations show a winter minimum in Rrs for the considered bands. However, different parameterizations of DOC-release induce a different accumulation of CDOM, especially in the eastern Mediterranean, and a different Rrs signature: a more active microbial loop during summer implies a decrease of Rrs at 412 nm. We demonstrate how the usage of a bio-optical model allows us to corroborate hypotheses on CDOM-cycling based on blue–violet Rrs data, supporting the importance of this complementary data stream with respect to satellite-derived chlorophyll.


Introduction
One of Earth's largest exchangeable reservoirs of carbon is the dissolved organic matter (DOM) in the ocean. The pool of dissolved organic carbon (DOC) approximates the amount of carbon resident in atmospheric carbon dioxide [1]. DOC plays a key role in the global carbon cycle and in the functioning of marine ecosystems [2]. A fraction of DOM is colored or chromophoric (CDOM), and absorbs light at the UV and visible wavelengths. It can also re-emit a fraction of absorbed light as fluorescence (FDOM). When CDOM absorption is high, water becomes yellow or brown [1]. CDOM optical properties (i.e., absorption and fluorescence) can thus help provide additional information about DOM pool dynamics and its main sources.
The concentrations of DOC in the Mediterranean Sea are similar to those observed in the open ocean, spanning between 31 and 128 µM, with values up to 366 µM in coastal areas impacted by river inputs [3]. DOC concentrations and CDOM absorption in surface waters can at times be higher than average values, causing a green hue of the water when observed from satellites [4]. Greener waters would generally imply higher chlorophyll concentrations than the ones measured in-situ [4,5]. As a consequence, satellite algorithms may provide higher chlorophyll concentration measurements than those taken in-situ for such greener waters. Such bio-optically anomalous behavior, in the Mediterranean Sea, has been ascribed to different causes: enhanced particle backscattering and non-algal absorption due to Saharan dust deposition [5], and/or the presence of coccolithophores [6,7]. Another possible cause could be an enhanced CDOM concentration at surface layers [8,9], resulting in twice the concentrations observed in the Atlantic Ocean [4].
The lack of direct measurements of CDOM stimulated several case studies to tackle CDOM spatio-temporal dynamics in the Mediterranean Sea. For example, higher CDOM absorption at 440 nm (and thus higher CDOM concentration) for a given Chl concentration was observed at two sites: at the BOUSSOLE site during oligotrophic periods (Chl less than 1 mg chla m −3 ) [10] and in the Catalan-Balearic Sea, based on data from summer cruises in 2011−2012 (for Chl between 0.04 and 0.1 mg chla m −3 ), described in [11]. Ref. [10] showed also that light absorption at 440 m in surface waters at the BOUSSOLE site was largely dominated by CDOM during all seasons except for the algal bloom in March and April. Ref. [12] suggested different surface and subsurface dynamics, with surface values correlating to sea surface temperature, whereas subsurface CDOM exhibited tighter correlations with chlorophyll concentration.
DOC and CDOM production can be originated as a byproduct of microbe metabolism in the ocean (autochthonous origin) or derive from river discharge and atmospheric deposition (allochthonous origin): in the present work we will focus on the autochthonous origin due to plankton excretion. Up to 90% of phytoplankton-dissolved exudates is made up of polysaccharides with smaller proportions of amino-sugars, acidic sugars, and chitin [13]. Although most of them are optically inactive, a fraction of the DOM pool absorbs light at both ultraviolet (UV) and visible wavelengths [14]. These exudates provide a highly diverse class of molecules that can be used by heterotrophic bacteria as substrates for metabolic processes [15] and CDOM can be also produced as a byproduct of DOM metabolism by bacteria [10]. Within these autochthonous sources, the relative contribution to CDOM production by phytoplankton and bacteria is still under debate [16].
The Mediterranean Sea present different productive regimes with regions exhibiting relatively high biomass accumulation [17] and Levantine ultraoligotrophic regions exhibiting strong phosphorus limitation [18,19]. Marginal Sea (i.e., the Adriatic Sea and Aegean Sea) are characterized by a strong influence of river discharge and in the case of the Aegean Sea by the interaction with Dardanelles. Microbial loop dynamics have been shown to be an important regulating mechanism for primary producer dynamics [20], as an alternative path with respect to the classical food chain. In particular, heterotrophic bacteria play an important role in terms of DOC remineralization, nutrient regeneration, and resource competition with autotrophic microbes [20]. The heterogeneity of the trophic gradients and the microbial loop dynamics can have an important role in determining DOC cycling and CDOM distribution in the Mediterranean Sea.
In this work we introduce remote sensing reflectance (Rrs) metrics to measure the skill of different CDOM formulations. Rrs is the ratio of the water-leaving radiance in the air to the downward plane irradiance incident onto the sea surface, with both measured just above the sea surface [21]. The importance of this quantity is related to the fact that it is the quantity measured by satellites and it is used, after atmospheric corrections, in computing all the ocean color products such as satellite chlorophyll. The model used in this work allows us to reconstruct such Rrs measurements based on the biogeochemical properties of the water mass considered and compare them directly with satellite data.
In this study we adopt a 3-dimensional biogeochemical model to investigate how different formulations of carbon cycling and CDOM production could affect chlorophyll and reflectance data and examine how such a modelling tool can be used to evaluate different physiological mechanisms and parameterizations. The manuscript is organized as follows: in the materials and methods section the biogeochemical model, the radiative transfer model, and the validation data are presented. Specific description is provided for each of the three alternative model configurations adopted. In the results section, model results are compared with satellite data in particular chlorophyll product and remote sensing data. Novel findings are examined in the discussion section and summarized in the conclusions section.

The Transport Reaction Model
The Mediterranean Sea biogeochemical model OGSTM-BFMv5 is based on the software detailed in [22,23] and is off-line coupled to the 1/16 • version of the Mediterranean Copernicus Marine Environmental Monitoring Service (CMEMS, http://marine. copernicus.eu) physical system [24].
The physical forcing data needed to compute the transport due to water mass circulation include the 3D horizontal and vertical current velocities, vertical eddy diffusivity, potential temperature, and salinity and 2D surface data for wind stress. These datasets were produced by the CMEMS Mediterranean Monitoring and Forecasting (MED-MFC) system.
The horizontal resolution is approximately 6 km and there are 72 vertical levels with 3 m resolution at the surface, reducing to 300 m at the deeper levels.
The transport module that accounts for advection and diffusion is coupled with a biogeochemical model that resolves the element cycling according to biogeochemical processes like, for example, carbon uptake by phytoplankton and mineralization by bacteria. Additionally, a multispectral radiative transfer model (RT model) resolves the propagation of light along the water column according to the water properties and the optically active constituents like chlorophyll, CDOM, and the particulate matter in it. Feedbacks are considered, so the biogeochemical properties affect light propagation that in turn affects biogeochemical properties modulating primary production and photo-acclimation of chloroplasts. All these processes are included, in a simplified form, in the model and are schematized in Figure 1. Schematics of the biogeochemical model interactions with the multispectral light propagation model. Phytoplankton, particulate, and chromophoric-dissolved organic matter (CDOM) (as a fraction of dissolved organic carbon (DOC)), are responsible for absorption and scattering (forward and backscattering). In turn the light affects the biogeochemistry modulating primary production and photobleaching. The model allows us to compute remote sensing reflectance (Rrs) at all the wavelengths considered (400 nm and 425 nm shown).

The Biogeochemical Model
The biogeochemical model here adopted (Biogeochemical Flux Model (BFM) [25]) has been used to study primary producers' biogeochemistry [17], alkalinity spatial and temporal variability [26], and CO2 fluxes [27], and has been corroborated using in-situ data for operational purposes within CMEMS [23]. A full description of the BFM model equations can be found in [25]. The BFM model has been expanded in the present configuration adding CDOM dynamics by assuming that in DOC production only a fraction of the produced DOC is optically active. The carbon dynamics in the BFM phytoplankton module

The Biogeochemical Model
The biogeochemical model here adopted (Biogeochemical Flux Model (BFM) [25]) has been used to study primary producers' biogeochemistry [17], alkalinity spatial and temporal variability [26], and CO 2 fluxes [27], and has been corroborated using in-situ data for operational purposes within CMEMS [23]. A full description of the BFM model equations can be found in [25]. The BFM model has been expanded in the present configuration adding CDOM dynamics by assuming that in DOC production only a fraction of the produced DOC is optically active. The carbon dynamics in the BFM phytoplankton module are determined by the rates of primary production, respiration, excretion, and predation. We assign for each production process of DOC the corresponding fraction of optically active molecules (CDOM): this is facilitated by the modular structure of BFM.
The photo-bleaching processes are based on the same formula introduced in [28]. Consistently with three pools of DOC (labile, semi-labile, and semi-refractory), the novel BFM model configuration includes three forms of CDOM, each CDOM component is degraded by heterotrophic bacteria according to its lability similar to the corresponding DOC variable. Moreover, the three forms of CDOM share the same response to light degradation with an e-folding time scale of one month.
Three numerical experiments were performed to evaluate the effect of the different parameterizations on model skill metrics shown in the results section. Here we report three parameterizations involving alternative CDOM and DOC cycling formulations ( Figure 2): • EXP-1: the CDOM is partitioned as 2% of DOC in all the fluxes [28]. The formulation includes a carbon correction formulation than in nutrient stressed conditions implies that the excess carbon synthesized is excreted in the form of semi-labile DOC [29,30]; • EXP-2: the CDOM is partitioned as 2% of DOC in all the fluxes but the nutrient stress limitation is in part accounted as non-realized production [31] and in part channeled to the production of semi-labile DOC. In this experiment we chose to remove the flux to semi-labile DOC to explore an extreme condition; • EXP-3: the CDOM production is related to activity exudation and is considered colored, the semi-refractory component of CDOM in not considered.

The Radiative Transfer Model
The light propagation is resolved coupling an atmospheric multispectral radiative transfer model (RT model) [32] with an in-water RT model [28,33], with bands of 25 nm resolution in the UV and visible wavelengths, as in Figure 1.
The equations considered in the model resolve three state variables: a direct down-

The Radiative Transfer Model
The light propagation is resolved coupling an atmospheric multispectral radiative transfer model (RT model) [32] with an in-water RT model [28,33], with bands of 25 nm resolution in the UV and visible wavelengths, as in Figure 1.
The equations considered in the model resolve three state variables: a direct downward component (Ed), a downward diffuse component (Es), and an upward diffuse component (Eu): ( The coefficients C, B, and F depend on the absorption and scattering of the biogeochemical constituents along the water column. In particular the light attenuation (in terms C) depends on water, chlorophyll concentration, particulate concentration, and CDOM concentration. The light absorption of CDOM is described using reference absorption at 450 nm of 0.015 m 2 (mgC) −1 [28] and an exponential slope of 0.017 nm −1 [10,34] (Figure 3). The three-stream equations correspond to a boundary value problem and the determination of the solutions requires information on Ed and Es just below the sea surface that are obtained using the results of the Ocean-Atmosphere Spectral Irradiance Model (OASIM, [32]). Radiometric data and absorption spectra of CDOM and phytoplankton can be used to constrain the bio-optical model and evaluate its skills. Absorption and scattering coefficients were adopted from the data set used in the NASA Ocean Biogeochemical Model (NOBM) coupled with the OASIM model provided by Dr. Watson Gregg (NASA). C) depends on water, chlorophyll concentration, particulate concentration, and CDOM concentration. The light absorption of CDOM is described using reference absorption at 450 nm of 0.015 m 2 (mgC) −1 [28] and an exponential slope of 0.017 nm −1 [10,34] ( Figure 3). The three-stream equations correspond to a boundary value problem and the determination of the solutions requires information on Ed and Es just below the sea surface that are obtained using the results of the Ocean-Atmosphere Spectral Irradiance Model (OASIM, [32]). Radiometric data and absorption spectra of CDOM and phytoplankton can be used to constrain the bio-optical model and evaluate its skills. Absorption and scattering coefficients were adopted from the data set used in the NASA Ocean Biogeochemical Model (NOBM) coupled with the OASIM model provided by Dr. Watson Gregg (NASA). Details concerning the RT model solution are provided as Supplementary Materials, the approach is founded on the analytical solutions presented in [28]. Simulated remote sensing reflectances (Rrs) are calculated following [28], i.e., from reflectance (R = Eu/(Ed+Es)) divided by a bidirectional factor Q [21]. Q is expressed in steradians and depends on time and latitude according to the equation introduced by [35].
Rrs is presented in two forms: as derived in [28] valid just below the sea surface, so called Rrs(0-), and computed considering a further correction [36,37] to derive Rrs(0+) that is valid above the sea surface. The latter is the most appropriate to be compared with satellite observations. In the following sections we keep both Rrs(0+) and Rrs(0-) to evidence the sensitivity of such correction in terms of Rrs determination. Details concerning the RT model solution are provided as Supplementary Materials, the approach is founded on the analytical solutions presented in [28]. Simulated remote sensing reflectances (Rrs) are calculated following [28], i.e., from reflectance (R = Eu/(Ed + Es)) divided by a bidirectional factor Q [21]. Q is expressed in steradians and depends on time and latitude according to the equation introduced by [35].

Dataset for Model Validation
Rrs is presented in two forms: as derived in [28] valid just below the sea surface, so called Rrs(0-), and computed considering a further correction [36,37] to derive Rrs(0+) that is valid above the sea surface. The latter is the most appropriate to be compared with satellite observations. In the following sections we keep both Rrs(0+) and Rrs(0-) to evidence the sensitivity of such correction in terms of Rrs determination.

Dataset for Model Validation
Remote sensing reflectance is the main satellite product for ocean color applications as it represents the primary variable for the derivation of other bio-optical and biogeochemical products, such as chlorophyll a concentration (total Chl and phytoplankton functional types (PFTs)), diffuse attenuation coefficient for downwelling irradiance (Kd490), as well as inherent optical properties (IOPs), e.g., phytoplankton light absorption aph(443), absorption of detrital and dissolved organic matter adg(443), and particulate optical backscattering bbp(443).
To extend the existing time series beyond that provided by a single satellite sensor, the European Space Agency (ESA) has in recent years promoted the Ocean Color-Climate Change Initiative (OC-CCI), a multi-sensor, global, ocean-color product, mainly devoted to climate applications [38], that merges observations from four different sensors: SeaWiFS, MODIS, MERIS, and VIIRS. Current operational algorithms are detailed in the CMEMS Quality Information Document [39], keeping SeaWiFS wavelengths as the reference ones.
The rationale for such a choice most likely stems from the fact that SeaWiFS is considered the highest quality sensor with the best match to in-situ observations up to date [40].
For the model validation we selected the following CMEMS products: •

Simulation Protocol
The simulations are based on the Med-BFM system used operationally to perform short term simulations in the Mediterranean Sea [23]. All the simulations presented have been spun up from 2004 to 2010 with the EXP-1 configuration and then the specific configurations of EXP-1, EXP-2, and EXP-3 have been launched from 2011 to 2017. The year 2017 is presented in the results, and is used to compute the statistics and metrics. For each experiment the 3D monthly average fields are analyzed to filter out the high frequency signals, which are not the focus of the present work. The statistical analysis of the model and satellite comparison is performed for the Mediterranean sub-basins of Figure 4. The subdivision, currently used to validate CMEMS operational products [23], reflects a synthesis of published bioregionalization of the Mediterranean Sea [17,42,43], and allows for within-region consistent average calculations. Thus, the spatial averages over each sub-basin are computed for the monthly model and satellite means, and Pearson correlation and root mean square difference (RMSD) are computed on the timeseries of the sub-basin model and satellite monthly means. State variables of all simulations have been initialized according to [17,19].
We also considered the riverine input of semi-refractory CDOM assigning a 2% ratio of organic carbon runoff (results not shown), but the effects on basin-wide scale were negligible given the short time scale of the simulation (seven years), and especially considering the long time scale of CDOM degradation adopted in the present simulations (i.e., a time scale of 90 years at 10 • C temperature). We noted that when using a mineralization time scale of semi-refractory CDOM of one year [28] the deep water inventory was unstable and was rapidly depleted. Besides river contribution, also atmospheric inputs of optically active substances where not considered, this implies that the results presented are basically driven by autochthonous CDOM. thesis of published bioregionalization of the Mediterranean Sea [17,42,43], and allows for within-region consistent average calculations. Thus, the spatial averages over each subbasin are computed for the monthly model and satellite means, and Pearson correlation and root mean square difference (RMSD) are computed on the timeseries of the sub-basin model and satellite monthly means. State variables of all simulations have been initialized according to [17,19]. We also considered the riverine input of semi-refractory CDOM assigning a 2% ratio of organic carbon runoff (results not shown), but the effects on basin-wide scale were negligible given the short time scale of the simulation (seven years), and especially considering the long time scale of CDOM degradation adopted in the present simulations (i.e., a time scale of 90 years at 10 °C temperature). We noted that when using a mineralization time scale of semi-refractory CDOM of one year [28] the deep water inventory was

Chlorophyll Spatial Distribution
Global and regional biogeochemical models are usually validated with respect to satellite ocean color (OC) products. The present model can provide chlorophyll both as model phytoplankton physiological units (i.e., state variable of the phytoplankton functional type) and/or as a derived quantity from Rrs output, by applying the same empirical formulas used by remote sensing. The model data comparison will be based on the model phytoplankton physiological unit [17] since it allows the direct investigation of phytoplankton dynamics and ecosystem functioning. While the model simulates specific chlorophyll for each of the four PFTs, we compared the total model chlorophyll with the corresponding satellite product. The dataset used for the comparison is the CMEMS satellite product that adopts a regional algorithm specifically calibrated for the Mediterranean Sea [39,44].
All the experiments (EXP-1, EXP-2, and EXP-3) produce a west-east negative gradient with higher chlorophyll values in the western basin consistently with the OC data ( Figure 5). The novel bio-optical setup produces a west to east deepening gradient of the deep chlorophyll maximum (DCM) depth, covarying with the surface chlorophyll negative gradient, with depths of approximately 30 m in the Alboran Sea sub-basin (alb, Figure 4) and 100 m in the Levantine basin (Figure 5c,e,g). In EXP-1 the chlorophyll concentrations at the DCM in the 10E-35E longitude part of the transect (Figure 5c) are in the range of in-situ measurements and appear slightly underestimated in the 5W-10E longitude part of the transect [17,45]. The EXP-2 simulation shows higher values of chlorophyll which can be related to the different microbial loop formulation. In fact, in this case, the overflow of carbon is fully directed to DIC rather than to DOC, impacting the substrate composition available for bacteria uptake. The reduced apportioning to DOC strongly depresses the microbial loop with a consequent advantage for phytoplankton, which shows more than double the accumulation especially at surface. Moreover, higher concentration values are found in the subsurface at the DCM depth. The chlorophyll distribution simulated in EXP-3 appears similar to EXP-1.
the overflow of carbon is fully directed to DIC rather than to DOC, impacting the substrate composition available for bacteria uptake. The reduced apportioning to DOC strongly depresses the microbial loop with a consequent advantage for phytoplankton, which shows more than double the accumulation especially at surface. Moreover, higher concentration values are found in the subsurface at the DCM depth. The chlorophyll distribution simulated in EXP-3 appears similar to EXP-1.

Remote Sensing Reflectance
This novel bio-optical model allows the reconstruction of remote sensing reflectance (Rrs) that determines the color of the water. Rrs can be visualized as a spectrum (Figure 6),

Remote Sensing Reflectance
This novel bio-optical model allows the reconstruction of remote sensing reflectance (Rrs) that determines the color of the water. Rrs can be visualized as a spectrum (Figure 6), that typically has higher values in the blue and green regions and declines towards the red bands.  Confirming what is generally known [46], from the extensive visual analyses of such Rrs spectra plots ( Figure 6) we understood that, in the open ocean, a large part of the variability of the reflectance appears in the shorter wavelengths considered (i.e., <500 nm) due to phytoplankton and CDOM optical contributions.
The <500 nm region shows the larger deviations between the model and satellite, therefore it makes sense to focus on this interval to define the metrics and to perform model results analyses.
In general, the water color simulated by the model is less blue compared to the Confirming what is generally known [46], from the extensive visual analyses of such Rrs spectra plots ( Figure 6) we understood that, in the open ocean, a large part of the variability of the reflectance appears in the shorter wavelengths considered (i.e., <500 nm) due to phytoplankton and CDOM optical contributions.
The <500 nm region shows the larger deviations between the model and satellite, therefore it makes sense to focus on this interval to define the metrics and to perform model results analyses.
In general, the water color simulated by the model is less blue compared to the CMEMS OC product; in fact, the model tends to underestimate the reflectance in the <500 nm interval with the exception of the marginal Seas, particularly for EXP-1 and less so for EXP-2 and EXP-3, the latter shown in Figure 6.
Marginal Seas (Adriatic and Aegean) and areas particularly influenced by coastal inputs (e.g., lev3 due to the Nile river) present different observed spectra with high reflectance variability in the range 500-600 nm. Among the sub-basins, also ion1 shows deviations in the 550 nm bin related to the shelf area of the Gulf of Gabes area where sea floor effects, not considered in the present radiative transfer model formulations, impact on the simulated Rrs.
In the marginal Seas, the satellite Rrs in the 500-600 nm range was higher than simulated pure water (not shown). This suggests that the higher reflectance is probably related to higher scattering by particles, probably of terrestrial origin. Therefore, since we did not consider such effects we preferred, for the present study, to focus on the <500 nm wavelengths. Clearly, the novel bio-optical model has the potential to investigate water masses with higher optical complexity, as the cases of the marginal Seas, however this would require a specific study that should consider additional factors such as terrestrial input of particulate substances, optically active DOM [47], and the effects of sea floor reflection.
Bio-optical constituent accumulation and consumption is strongly correlated with seasonal cycling, an example is the excretion of CDOM modulated by primary production. Therefore, information related to average Rrs spectra, as shown in Figure 6, can be augmented by analyzing its seasonally variability. In particular, considering that much of the bio-optical absorption takes place in the blue wavelengths, the 412 nm band is the most informative radiometric parameter available in the considered validation dataset.
During winter months (January, February, and March) all the tests show a minimum in Rrs (Figures 7-9). In the western sub-basins (alb, swm1, swm2, nwm, tyr1, and tyr2) the minimum of Rrs co-varies with the maximum in chlorophyll and CDOM (Figure 7), both for the EXP-1 experiment and satellite data. In the same sub-basins, the Rrs maximum is measured by satellite in the summer while it peaks in October/November in the EXP-1 model corresponding to the minimum in CDOM surface concentration. In the eastern sub-basins (ion1, ion2, ion3, lev1, lev2, lev3, and lev4), an Rrs minimum is also observed in winter (with one month simulated delay) but strikingly EXP-1 simulated a second minimum during the summer/early fall that is absent in the observed radiometry data. This minimum in reflectance co-varies with CDOM dynamics: it is mirrored by a maximum of CDOM in the same period. EXP-1 simulated an accumulation pattern of CDOM in the eastern sub-regions and especially in the Levantine area due to ultraoligotrophic nutrient limitation that produced a minimum of Rrs in late summer/start of fall.
basins (ion1, ion2, ion3, lev1, lev2, lev3, and lev4), an Rrs minimum is also observed in winter (with one month simulated delay) but strikingly EXP-1 simulated a second minimum during the summer/early fall that is absent in the observed radiometry data. This minimum in reflectance co-varies with CDOM dynamics: it is mirrored by a maximum of CDOM in the same period. EXP-1 simulated an accumulation pattern of CDOM in the eastern sub-regions and especially in the Levantine area due to ultraoligotrophic nutrient limitation that produced a minimum of Rrs in late summer/start of fall.  [41], thick blue lines (model results from EXP-1 simulation just above the sea surface (0+)) and thick cyan lines (EXP-1 simulation just below the sea surface (0-)). Data and model are averaged over each sub-basin (see Figure 4 for definition). Seasonal variability of the total surface chlorophyll (units mgChl/m3): green crosses (CMEMS product [39]) and green lines (model output) scale is on the left y-axis. Black dots represent the modelled seasonal evolution of total CDOM concentration normalized to one.     In the EXP-2 simulation, the winter minima in Rrs are simulated in all the sub-basins while a peak of reflectance is simulated in summer. The reflectance is generally underestimated in the whole experiment. Contrary to the EXP-1 simulation, the chlorophyll of EXP-2 appears to be the controlling variable driving Rrs. The winter minimum in reflectance starts when chlorophyll accumulates in the fall, this accumulation is much higher than the chlorophyll satellite estimate (Figure 8). In EXP-2 the reduced content of CDOM due to reduced microbial loop components has a clear impact both on chlorophyll and reflectance. The increase of chlorophyll does not compensate the reduction of CDOM in terms of absorption and the reflectance is less underestimated in all the sub-basins than in EXP-1.
The modulation of microbial loop cycling in term of Rrs and chlorophyll is non-linear, as is expected given the intrinsic non-linear nature of biogeochemistry. In fact, considering an intermediate configuration between EXP-1 and EXP-2 where half of the carbon overflow due to nutrient stress is channeled in the microbial loop as DOC, and the other part to DIC, the model simulates chlorophyll variability similar to EXP-1 but Rrs similar to EXP-2 (results not shown).
The differences in the model behavior both for chlorophyll and reflectance observed in EXP-1 and EXP-2 are related to the different parameterizations of DOC exudation considered, this in turn is reflected in the CDOM production.
The two examples shown in Figures 7 and 8 identify alternative regimes of resource competition between phytoplankton and bacteria. In EXP-1 the polysaccharides exudate released by primary producers due to nutrient-stressed conditions allows bacteria to compete for the low nutrient availability. In EXP2, the reduced exudation parameter channels less carbon in the microbial loop therefore damping bacteria and in turn favoring phytoplankton biomass.
The fact that such increased phytoplankton absorption provides the best magnitude and seasonality of Rrs indicates that whatever it is shaping the Rrs, it should co-variate with phytoplankton net growth. This was the reasoning we followed to define EXP-3, where the release of CDOM is synchronized with activity and net growth of phytoplankton.
The EXP-3 shows a correct chlorophyll seasonality with only an overestimation in winter months but more importantly the Rrs has a correct seasonality (Figure 9), with a minimum of reflectance in winter and a peak of reflectance in summer. All the sub-basins show a reduction in Rrs in late fall with the exception of lev1, lev2, lev3, and lev4. In the latter, the summer maxima of reflectance persist till late fall and early winter. The fall reduction of Rrs co-varies with an increase in CDOM, in turn synchronized with the increase in chlorophyll.
The summary statistics of the three experiments (Table 1) confirm the results derived from the plot analysis. In particular, EXP-1 has no significant correlation for all the three wavelengths in all the sub-basins with the exception of the western ones (alb, swm2, and nwm), and it has relatively low RMSD on simulated chlorophyll. In the EXP-2 experiment the correlation is significant in all the sub-basins with high correlation demonstrating that the modelled seasonality is coherent to that inferred from space. The only exceptions are the northern Adriatic Sea (adr1) for Rrs at 412 nm and 443 nm, and the southern Adriatic Sea (adr2) for Rrs 490. As noted before, marginal Seas are characterized by optically complex waters and require specific analysis. In EXP-2 the RMSD of chlorophyll is much higher than in the case of EXP-1. Interestingly, marginal Seas do not appear particularly faulty as, on the contrary, can be noted by examining the Rrs for the same experiment. EXP-3 appears to be a good compromise between EXP-1 and EXP-2: in fact the correlation of Rrs is significant almost everywhere with the exception of the marginal Seas and the RMSD of chlorophyll is comparable to EXP-1.
The temporal covariance between Rrs and chlorophyll is consistent with the approaches usually adopted in empirical models of ocean color. In particular, in these algorithms the accumulation of chlorophyll is determined by a reduction of Rrs (at a given wavelength according to the specific algorithm) normalized versus a reference measured Rrs at 555 nm (a sort of blank) that is less affected by chlorophyll absorption. In our results we focus on Rrs at 412 nm because it is most interesting for CDOM dynamics and we observe similar dynamics: an increase in phytoplankton and CDOM results in an increase in absorption and a corresponding reduction in Rrs. Table 1. Validation metrics of Rrs and chlorophyll for the three experiments. Data and model are spatially aggregated for the year 2017 over the sub-basins shown in Figure 4. Metrics are temporal correlation for Rrs(0+) at 412 nm, 443 nm, and 490 nm, with respect to CMEMS Rrs ocean color (OC) product [41], and RMSD (units mgChl/m 3 ) with respect to CMEMS satellite chlorophyll product [39]. Entries in bold indicate that correlation is significant with p-value < 0.05.

Discussion
The experiments illustrated in the previous section show how the usage of the biooptical model offers the opportunity to validate different CDOM formulations according to radiometric-based metrics.
Both EXP-1 and EXP-2 present some flaws. With the standard metrics based on chlorophyll [23], EXP-1 clearly appears the one with lower RMSD. On the contrary, in terms of Rrs at 412 nm, EXP-2 shows a better skill even if in this case the chlorophyll is overestimated and the microbial loop is reduced too much. As evidenced by EXP-1 and EXP-2, the production of DOC, and in turn CDOM, is not synchronized with net primary production and with phytoplankton biomass accumulation. During the winter period, phytoplankton is the mayor driver of the system and competition with bacteria is less important, as demonstrated by similar simulated winter chlorophyll blooms between EXP-1 and EXP-2. During summer the availability of DOC, in particular the semi-labile component, activates bacteria growth and enhances the uptake of bacteria of the limiting nutrients that strongly limit phytoplankton growth. The important role of the microbial loop in the ultraoligotrophic areas of the Mediterranean Sea has been evidenced by experimental and theoretical studies [20].
Moreover, similar results to the ones here presented were obtained by zero-dimensional modeling experiments [48] and a 3D modelling study for the Adriatic Sea [49] where two regimes were identified: the winter regime when phytoplankton dominates production and a summer regime when the bacteria are able to outcompete phytoplankton due to large release of DOC. It was demonstrated how the excretion of large pools of semi-labile DOC triggers the bacterial uptake of available nutrients and limits, through resource competition, phytoplankton growth.
Model results indicate that to reproduce realistic Rrs data the bacteria metabolization of excreted DOC due to nutrient limitation should not be a major source of CDOM, or at least the fraction CDOM:DOC should be much smaller than the 2% that is used in the case of EXP-1 and EXP-2.
The only sub-basin showing a marked minimum of reflectance in summer is the northern Adriatic Sea, where the high productivity triggered by riverine discharge [49] could strongly affect the optical properties of the water and possibly substantially trigger bacterial CDOM production in summer. The northern Adriatic Sea is highly influenced by terrestrial inputs and it is characterized by optically complex waters, therefore additional analyses will be needed to evaluate model results.
Among the experiments performed, EXP-3 appears to be the best configuration with good skill in both the indicators (i.e., Rrs correlation and chlorophyll RMSD). EXP-3 differs in two aspects with respect to EXP-1 and EXP-2: the release of CDOM in the form of labile component is synchronized to the activity excretion of phytoplankton and to bacteria mortality, and the 2% rate of semi-refractory CDOM:DOC ratio produced by heterotrophic bacteria is not included. EXP-2 shows that the synchronization of absorption with phytoplankton cycling improves reflectance and we assume that nutrient replete conditions (during production period) are likely necessary to produce optically active substances.
The reduced production of semi-refractory CDOM (EXP-2 and EXP-3) always produces an underestimation of reflectance. Clearly the parameterization of semi-refractory CDOM components requires further in-depth analyses considering in-situ measurements (e.g., [50]).
In all the experiments performed, the deep ocean inventory of CDOM is kept invariant. This indicates that in terms of seasonal variability of Rrs it is less relevant as compared to the differences among EXP1-3. Even if we focused on surface dynamics, the deep inventory of CDOM could be relevant due to mixing or upwelling. Clearly, long term simulations should be performed to reconstruct the dynamic of the refractory component of CDOM.
The formula to derive Rrs(0+) from Rrs(0-), accounting for water/air interface effects on irradiance, applied to model results based on the formula by [36], strongly impacted the absolute value of Rrs switching from the overestimation to underestimation of corresponding remotely sensed reflectance. Given the low value of reflectance, the correction is roughly equivalent to dividing Rrs by approximately a factor of two. This correction does not impact the seasonal variability of the signal but changes the model solutions from a condition of overestimation in the blue to an underestimation of Rrs. Moreover, it introduces an underestimation in the Rrs at 490 nm. The formula is based on solutions determined with the state-of-the-art Hydrolight software [51], and it would important to evaluate if specific corrections are needed for an optically anomalous area like the Mediterranean Sea and/or if the discrepancy in the Rrs results can be related to intrinsic approximations of the three-stream approach.
An interesting aspect of the present simulations is that the explicit definition of different physiological processes in the Med-BFM biogeochemical model allows us to assign them specific optical properties (CDOM fraction) and to evaluate their impact on a measurable quantity such as Rrs. Clearly during summer a significant part of exudate is transparent: in fact the high level of DOC production is not associated with a minimum of Rrs, as evident from satellite data. This is consistent with several experimental results where it is shown that under carbon overflow the exudate will be dominated by transparent carbohydrates [52].
The fact that DOC and CDOM formation follows different pathways clearly indicates that it is, in principle, difficult to extrapolate and correlate one in terms of the other only on the basis of the simultaneous concentration distributions. From a study conducted in the Sargasso Sea [53] it emerged that CDOM is not a direct product of primary producers given the spatial and temporal decoupling between chlorophyll and CDOM observed. On the contrary, from the analyses reported here it seems that for the Mediterranean Sea there is a certain synchronization between primary producers' cycles and CDOM formation, at least for the surface waters, with peaks in winter of chlorophyll and CDOM and corresponding minima in Rrs at 412 nm.
A time shift in the north western Mediterranean Sea, mediated by microbial digestion of phytoplankton with approximately one month lag between chlorophyll peak and CDOM maxima, was documented [10]. To analyze in more detail the sub-basin dynamics, a high frequency modelling study could be performed to better investigate synchronization.
Synchronization between chlorophyll and CDOM production has been documented in the North Atlantic subtropical gyre [30], where characteristics of the area could present similarities to the oligotrophic regions of the Mediterranean Sea.
On the other hand, the role of heterotrophic bacteria in shaping CDOM spatial and temporal gradients has not emerged clearly in the results presented. Some indication of CDOM increase and a reduction of Rrs are found in model results in the eastern sub-basins in correspondence of summer bacterial activity, however, these findings are not strikingly confirmed looking at corresponding satellite Rrs data (Figure 9).
In the present work, we mainly focused on model versus satellite products, but increasingly abundant data are available from robotic autonomous profilers (BGC-Argo floats) [54], that measure physical and biogeochemical properties including bio-optical parameters. These instruments can be very useful for data driven models [55] and for validation [23,32]. The advantage of these tools is that they provide information along the water column and this information could be useful to understand CDOM dynamics below the sea surface [9,30]. BGC-Argo floats are equipped with sensors measuring planar downward irradiance at 412 and 490 nm, and in particular at 380 nm where CDOM absorption is sharply dominant with respect to phytoplankton. This additional information could be important to expand the scope of the present study and will be addressed in future works.

Conclusions
CDOM is an important biogeochemical property that strongly affects the bio-optical characteristics of the Mediterranean Sea. We show that with a medium complexity bio-geochemical model describing CDOM dynamics coupled with a multispectral radiative transfer model it is possible to analyze and discuss hypothesis on CDOM production and consumption. Different physiological formulations impact in a different way simulated chlorophyll and remote sensing reflectance allowing us to determine which formulation is the most consistent with observations. Moreover, the use of remote sensing reflectance appears extremely useful to discriminate marginal Seas and coastal areas versus open ocean areas given the remarkable differences in the spectral signature. All experiments show similar winter ocean color signals, with a minimum in Rrs for the considered bands. The most interesting differences are found in summer according to the different mechanisms of DOC and CDOM pathways chosen for each experiment. EXP-3, characterized by CDOM formation produced by phytoplankton activity excretion and no flux of CDOM to the semi-refractory pool, appears the most consistent considering chlorophyll and remote sensing reflectance data.