Diel variability affects the inorganic carbon system in the sea-surface microlayer and influences air-sea CO₂ flux estimates

The ocean plays a crucial role in the global carbon cycle by absorbing and storing about one-third of anthropogenic carbon dioxide (CO₂). It is estimated that the ocean has sequestered approximately 26 % of CO₂ emissions over the last decade, resulting in significant changes in the marine carbon system and impacting the marine environment. The sea-surface microlayer (SML) plays a crucial role in these processes, facilitating the transfer of matter and energy between the ocean and the atmosphere. However, most studies on the carbon cycle in the SML have primarily addressed daily variability and overlooked nocturnal processes, which may lead to inaccurate global carbon estimates. We analysed temperature, salinity, pH_T25, and pCO₂ using data collected over three complete diel cycles during an oceanographic campaign along the Croatian coast near Šibenik in the Middle Adriatic. Our analysis revealed statistically significant differences (p<0.05) between daytime and nighttime measurements of temperature, salinity, and pH_T25. Diel differences in pCO₂, were also observed, with patterns largely driven by temperature effects and short-term mixing. These differences may be related to the occurrence of buoyancy fluxes, which are typically more pronounced during the day and could enhance CO₂ fluxes, as observed with values of 1.98 ± 2.52 mmol cm⁻² h⁻¹ during the day, while at night, they dropped to 0.01 ± 0.02 mmol cm⁻² h⁻¹. These findings emphasise the importance of considering complete diurnal cycles to accurately capture the variability in thermohaline features and carbon exchange processes, thereby improving our understanding of the ocean's role in climate change.

1 Introduction

The ocean is a crucial climate regulator that mitigates the effects of anthropogenic emissions (Gattuso et al., 2015). It absorbs more than 90 % of the Earth's excess heat (Hoegh-Guldberg et al., 2014; Pörtner et al., 2019) and captures approximately one-third of anthropogenic carbon dioxide (CO₂) emissions (Pörtner et al., 2019; Wong et al., 2014). It is well known that CO₂ has a significant impact on seawater chemistry (Doney et al., 2009; Gattuso et al., 2015). It is estimated that the ocean has sequestered approximately 26 % of global CO₂ emissions over the last decade (Friedlingstein et al., 2025). Consequently, changes in the ocean environment have exceeded the magnitude and rate of natural variation due to anthropogenic carbon perturbation over the last millennia (Gattuso et al., 2015). In this context, the marine carbon system undergoes substantial modifications, resulting in various effects on the marine environment, including a decrease in pH levels (Doney et al., 2009; IPCC, 2023; Orr et al., 2005). Therefore, understanding the evolving state of marine biogeochemistry in the context of climate change is crucial, with a primary focus on the upper layers of the water column, which are closely linked to ocean-atmosphere interactions.

In the ocean-atmosphere system, the sea-surface microlayer (SML) plays a vital role in transferring materials and energy, such as heat, gases, and particles (Wurl et al., 2019), which must pass through it (Frka et al., 2009; Stolle et al., 2020). The SML is a distinctive and complex marine environment that represents the interfacial boundary layer between the ocean and the atmosphere (Cunliffe et al., 2013; Stolle et al., 2020; Wurl et al., 2011). Its thickness typically does not exceed 1000 µm, and it exhibits distinct biological and physicochemical properties compared with the underlying water masses (Cunliffe et al., 2013; Stolle et al., 2020; Wurl et al., 2011). Thus, the SML experiences instantaneous meteorological forcing, such as solar radiation, wind, and atmospheric inputs (Gassen et al., 2023; Wurl et al., 2019), which impacts the development of physical and biogeochemical processes occurring in the underlying water (Engel et al., 2017; Stolle et al., 2020; Wurl et al., 2017). These characteristics of the SML overlap with the growing interest in oceanographic research to study the spatiotemporal variability of the marine carbon system (Cantoni et al., 2016), as this layer is essential for understanding global marine biogeochemistry. However, to gain a comprehensive spatiotemporal perspective on marine carbon chemistry processes, focusing on underexplored fields, such as the role of nocturnal processes within the diel cycle, is crucial.

Daily variations force cyclic changes in chemistry (De Montety et al., 2011), which are influenced by processes such as photosynthesis, air-sea gas exchange, and various environmental conditions (e.g., light, temperature, and nutrient availability) (Poulson and Sullivan, 2010). These processes directly affect the seawater pH by adding or removing CO₂ from seawater (Poulson and Sullivan, 2010; Takahashi et al., 2002). It is well known that photosynthesis consumes CO₂ during the day, thereby reducing pCO₂ levels and consequently causing an increase in pH (Cantoni et al., 2012; Takahashi et al., 2002). At night, the CO₂ produced by respiration tends to be more constant and accumulates in the water column, leading to an increase in pCO₂ and a decrease in pH (Cantoni et al., 2012; del Giorgio and Williams, 2005; Gattuso et al., 1999; Shaw et al., 2012). Although these processes are significant for seawater chemistry, research has predominantly focused on diurnal processes. As a result, the roles of respiration and other nocturnal processes in the SML remain largely unexplored (Yates et al., 2007).

Within this framework, the Mediterranean Sea presents a unique research opportunity to investigate the spatiotemporal variability of the marine carbon system. It is often referred to as a “laboratory basin” (Bergamasco and Malanotte-Rizzoli, 2010; Robinson and Golnaraghi, 1994) because it enables us to approximate processes occurring on a global scale within a shorter timeframe and smaller space (Álvarez et al., 2014). Despite representing 0.8 % of the global ocean surface (Álvarez-Rodríguez, 2012), the Mediterranean Sea is considered an important anthropogenic carbon storage (Álvarez et al., 2014), as it absorbs a disproportionally larger amount of anthropogenic carbon compared to the global ocean (Hassoun et al., 2015; Schneider et al., 2010). Higher acidification ranges (−0.001 to −0.009 pH unit yr⁻¹) (Hassoun et al., 2022) were measured and exceeded those measured in the Atlantic Ocean (−0.001 to −0.0026 pH unit yr⁻¹) (Takahashi et al., 2014), in line with regional coastal studies showing strong acidification trends in the Mediterranean (Kapsenberg et al., 2017). However, as explored in this study, the dynamics of biogeochemical processes in coastal regions are more complex than those in the open ocean (Borges, 2005). When examining oceanic regions, it is essential to recognise that daily variability has a profound impact on marine carbon chemistry.

This study emphasises the role of diel variability in thermohaline features and dynamics of the marine inorganic carbon cycle in the coastal region of the Mediterranean Sea (Šibenik, Croatia), thereby providing a comprehensive understanding of the processes that influence the complete diurnal cycle. Therefore, it is essential to incorporate nocturnal processes into global estimates of marine carbon cycle dynamics. We present high-resolution data to understand biogeochemical processes and their variability in the water column, which complicates predictions of variations in the coastal marine carbon system (Cantoni et al., 2012). We have focused on understanding the nocturnal processes that affect the inorganic carbon system, thereby helping to clarify the uncertainties associated with this system. This study contributes to the identification of anthropogenic influences on the marine environment in the context of climate change.

2 Materials and Methods

2.1 Sampling Strategy and Seawater Analyses

Sampling was conducted in the lower estuary of the Krka River (central-eastern Adriatic), offshore in the St. Anthony Channel, near Šibenik (Fig. 1). The estuary is highly stratified and exhibits microtidality, extending approximately 23 km inland, with depths increasing from less than 2 m at the head to around 42 m at the mouth (Cukrov et al., 2024a; Prohić and Juračić, 1989). Similar to other Mediterranean estuaries, tidal ranges are small (0.2–0.5 m), resulting in minimal currents and stable vertical stratification, with freshwater or brackish water flowing over a lower marine layer (Cukrov et al., 2024a). Average annual river discharge is ∼ 50 m³ s⁻¹, varying seasonally from ∼ 5 to 480 m s⁻¹ (Bužančić et al., 2016; Cukrov et al., 2024b; Marcinek et al., 2020). Previous studies in the area have reported that water residence times range from a few days in winter to several weeks in summer, depending on the hydrodynamics of the freshwater or marine layers (Cetinić et al., 2006; Zutic and Legovic, 1987). During the sampling period (10–15 August 2020), precipitation was very low, averaging 1.4 mm d⁻¹ during the days preceding the campaign (Croatian Meteorological and Hydrological Service, 2020), consistent with the dry summer conditions typical for the region. This environment provided ideal stable conditions for assessing the daily variability of SML and ULW parameters with minimal interference from tides or precipitation.

Figure 1 Map of the sampling area in the Middle Adriatic Sea, created using Ocean Data View (Schlitzer, 2022). The red dot indicates the sampling area, whereas Zlarin Island, where the meteorological station is situated, is outlined in red.

Temperature, conductivity, and pH data from the SML were collected using sensors integrated in a flow-through system on the “Sea Surface Scanner (S³)” (Ribas-Ribas et al., 2017). The S³ is a 4.5 m long and 2.2 m wide uncrewed catamaran, remotely piloted by radio control from a small support vessel (a Zodiac). As specified in Ribas-Ribas et al. (2017), the system features rotating glass discs on which the SML water adheres via surface tension. A set of scraping mechanisms on the immersed side collects the adhering water and transfers it to a closed, continuous-flow system, minimising atmospheric contact. Underlaying water (ULW) was continuously collected at a depth of 1 m through a rigid inlet pipe. Both SML and ULW water were pumped directly to two independent packages of onboard sensors (Table 1) and sampler systems, each equipped with its own temperature, conductivity/salinity, and pH sensors, allowing simultaneous measurements to be taken in both layers. In addition, the S³ carried an automated water collector (Table 1) with twenty-four 1 L polypropylene bottles, which could be filled remotely with SML or ULW water, allowing discrete samples to be collected during transects. Discrete samples for inorganic carbon parameters (DIC and TA) were collected in bottles (250 mL) directly from the S³ outlets at each sampling event, with one sample per parameter and depth, approached by the catamaran with the support vessel, and immediately sealed to prevent atmospheric contamination.

Table 1 Manufacturers, models, and specifications of the sensors employed to measure pH, salinity, and temperature in the sea-surface microlayer (SML) and underlying water (ULW). Adapted from Ribas-Ribas et al. (2017).

Sensors recorded temperature, conductivity, and pH from both SML and ULW at 30 s intervals. Specifically, the pH was measured using the integrated sensor of the S³ flow-through system (Ribas-Ribas et al., 2017) and is reported on the total scale. The uncertainties of the measured temperature, conductivity and pH were estimated based on sensor specifications (Table 1). Solar radiance and wind data were collected from a meteorological station (Davis Instruments, Vantage Pro2 Plus) located on Zlarin Island. Furthermore, subsamples from the SML and ULW water collected by the S³ system were stored in high-density polyethylene (HDPE) bottles to determine the phosphate (PO${}_{\mathrm{4}}^{\mathrm{3}-}$) and silicate [Si(OH)₄] concentrations. To preserve the samples, mercury chloride (HgCl₂) was added. The preserved samples were stored at +4 °C until further analysis in the laboratory. Nutrient concentrations were measured using a sequential automatic analyser (SAA, SYSTEA EASYCHEM) according to the standard protocols described by Laskov et al. (2007) and Fanning and Pilson (1973). These nutrient data were subsequently used as inputs for calculating partial pressure of CO₂ (pCO₂) with the CO2Sys program (Version v3.2.0, MATLAB) (Sharp et al., 2020; Van Heuven et al., 2011).

Field deployments of the S³ were conducted six times per diel cycle between 10 and 15 August 2020, with each deployment lasting 30–45 min. Although the catamaran moved along short transects within the sampling area, for analytical purposes, the sampling location was treated as a single point. The catamaran remained within a small area (maximum horizontal displacement ∼ 90 m based on GPS positions), and consequently, no significant horizontal gradients were expected at this spatial scale, since spatial variations within the region were minimal. For the analysis, a total of three diurnal cycles were sampled. In each cycle, four deployments were conducted during daytime conditions (diurnal data), and two deployments were performed at night (nocturnal data). Environmental conditions during the campaign were typical for the region, with no significant rainfall events. This deployment strategy enabled representative sampling of the study area while minimising the potential influence of microtidal variations (Fig. S1 in the Supplement) on the SML and ULW measurements.

2.2 Marine Carbon System Determination

DIC samples (20 mL) were analysed by coulometric titration (CM 5014, UIC, USA) with an excess of 10 % phosphoric acid. Total Alkalinity (TA) samples (100 mL) were measured via potentiometric titration (916 Ti-Touch, Metrohm, Switzerland) and calculated using a modified Gran plot approach implemented in Calkulate (version 1.8.0) (Humphreys et al., 2022). Calibration was performed with certified reference material (Batch 187) obtained from A. G. Dickson at the Scripps Institution of Oceanography. The 1σ measurement precision was ±3 µmol kg⁻¹ for DIC and ±2 µmol kg⁻¹ for TA.

To eliminate the effects of temperature variation on the results, the pH values measured in the seawater samples were adjusted to the average temperature value of 25.43 °C, hereafter referred to as pH_T25. The normalisation utilised an adjustment factor of 0.018 pH units per °C, which is widely accepted for the range of pH and TA conditions typical of seawater (Eq. 1) (Dickson et al., 2007; Dickson and Millero, 1987; Zeebe and Wolf-Gladrow, 2001):

$$\begin{array}{}\text{(1)}& {\mathrm{pH}}_{T\mathrm{25}}=\mathrm{pH}+\left(T-\mathrm{25.43}\right)\cdot \mathrm{0.018}\end{array}$$

To evaluate the marine carbon system, pCO₂ was calculated with the CO2SYS program (Version v3.2.0, MATLAB) (Sharp et al., 2020; Van Heuven et al., 2011) using as input parameters: DIC; TA, PO${}_{\mathrm{4}}^{\mathrm{3}-}$, and Si(OH)₄; salinity (S); temperature (T), and pressure. The respective dissociation constants were used for carbon (Mehrbach et al., 1973), sulfate (KSO₄) (Dickson and Millero, 1987), fluorine (KF) (Perez and Fraga, 1987), and the borate-salinity ratio (Lee et al., 2010). Missing pCO₂ values were due to erroneous DIC and/or TA measurements. Consequently, standard deviations could not be calculated when the number of reliable data points was less than three. Once the pCO₂ values were calculated, the CO₂ flux through the ocean-atmosphere interface in the monitoring area was estimated using the following equation:

$$\begin{array}{}\text{(2)}& F=\mathrm{\Delta }p{\mathrm{CO}}_{\mathrm{2}}\cdot k\cdot \mathit{\alpha }\end{array}$$

where F is the CO₂ flux (mmol m⁻² d⁻¹), ΔpCO₂ is equal to the difference in the partial pressures of the gas between the surface water and the atmosphere, k is the gas transfer coefficient (cm h⁻¹) from Wanninkhof, (2014), and α is the solubility of CO₂ in seawater (mol L⁻¹ atm⁻¹). To calculate the partial pressure gradient of CO₂, atmospheric pCO₂ data were obtained as the monthly mean for August 2020 from the ICOS Lampedusa atmospheric station, which provides quality-controlled dry-air CO₂ measurements representative of the Mediterranean region (Pecci et al., 2023). To quantify the bias introduced when nocturnal fluxes were ignored, we compared daily means including all valid data from each cycle with means based only on daytime measurements. The percentage error was defined as $\left(\frac{F_{\mathrm{daytime}}-F_{\mathrm{diel}}}{F_{\mathrm{diel}}}\right)\times \mathrm{100}$, with both means calculated as arithmetic averages of available fluxes. This approach explicitly evaluates the uncertainty in daily estimates when nocturnal processes are excluded (Garbe et al., 2014).

2.3 Salinity and density correction

To ensure high-quality data, a correction factor (CF) was applied to the continuous salinity and pH_T25 data. Discrete salinity values obtained from laboratory analyses served as reference points and were compared to the average continuous values recorded by the S³ (Ribas-Ribas et al., 2017). This process resulted in the derivation of distinct CFs, which were applied at each depth and during each time interval of the S³ measurements (Ribas-Ribas et al., 2017). Once the salinity values were corrected, the pH was calculated using the CO2Sys program (Version v3.2.0, MATLAB) (Sharp et al., 2020; Van Heuven et al., 2011). These calculated pH_T25 values were then utilised as reference points for comparison with the average continuous values obtained from the S³ measurements (Ribas-Ribas et al., 2017). After correcting for salinity, the density (ρ) was calculated using the TEOS-10 (https://www.teos-10.org/index.htm, last access: 12 February 2025) equation of state in RStudio (RStudio Team, 2023) based on the observed temperature and salinity values. From this calculation, sigma-t was defined as the density at a given temperature and salinity minus 1000 kg m⁻³. The propagated uncertainties for sigma-t and pCO₂ were constrained by the measurement precisions of temperature (±0.1 °C), salinity (±0.2 %), pH (±0.005), together with the analytical reproducibility of DIC, and TA. These estimates provide a quantitative basis for evaluating the reliability of calculated variables in subsequent analyses.

2.4 Evaporation rate calculation

The evaporation rate (E) was estimated using the following formula, which relates the latent heat flux (Q_E) (Brutsaert, 2013), the latent heat of vaporisation (Lv) (Kittel and Kroemer, 1980), and the calculated density of seawater (ρ):

$$\begin{array}{}\text{(3)}& E={\displaystyle \frac{Q_E}{\mathrm{Lv}\cdot \mathit{\rho }}}\end{array}$$

2.5 Statistical Analysis

Since the assumptions of normality and homoscedasticity were not met for the collected data, the Kruskal–Wallis test was chosen as a nonparametric alternative to ANOVA. This test revealed significant differences between the data collected during the day and night, as well as between the SML and ULW, and among the medians of the cycles studied. A significance level (α) of 0.05 was established to determine whether the groups differed significantly. All statistical analyses were performed using RStudio (RStudio Team, 2023). In addition, to carry out a complete analysis of the variability of the marine carbon system during the study period, anomalies in temperature, salinity, pH_T25, and pCO₂ data were calculated using the following expression: Δ(SML – ULW). These differences were calculated for diurnal and nocturnal data for each cycle.

3 Results

The primary objective of this study was to understand the typical patterns of diurnal variability and depth-related differences in marine biogeochemistry during the observed cycles. Data on temperature, salinity, pH_T25, and pCO₂ were used to perform a statistical analysis of the similarity between the data collected during the diel cycle and to calculate the differences between the values measured in the SML and ULW. Additionally, we examined the temporal distribution using box-and-whisker plots and calculated the air-sea CO₂ exchange across the SML.

3.1 Meteorological conditions

To study the potential variance in the meteorological forcing observed during the observations, time-series graphs were plotted for solar radiation and wind speed (Fig. 2). A consistent pattern of solar radiation was observed during all three cycles, with peaks of 525, 404, and 420 W m⁻² observed at 14:00 UTC. Relatively low wind speeds were recorded during the day in all three cycles, averaging 1.26 ± 1.46 m s⁻¹. Maximum wind speeds were observed at 14:00 and 18:00 UTC, higher for Cycles 1 and 2 (2.14 and 5.22 m s⁻¹) than during Cycle 3 (0.95 m s⁻¹). Additionally, a decrease in wind speed was observed throughout the night, approaching near-zero values, except during Cycle 3, when the wind speed remained consistently close to zero.

Figure 2 Time series showing wind speed and solar radiation on Zlarin Island (Fig. 1) during the three studied cycles.

3.2 Daily trends

To investigate the overall distribution of physicochemical parameters throughout the daily sampling period at the SML and ULW, we show box-and-whisker plots for temperature, salinity, and pH_T25 across different time intervals (Fig. 3). We also performed the Kruskal-Wallis test to determine statistically significant differences (see Sect. 2.5, Materials and Methods) between the SML and ULW (Table S1). The analysis revealed no significant temperature differences between the two depths during Cycle 1 (SML: 25.00 ± 0.85 °C; ULW: 24.90 ± 0.72 °C) and Cycle 2 (SML: 25.40 ± 0.62 °C; ULW: 25.40 ± 0.83 °C). However, in Cycle 3, ULW was slightly warmer, with a mean temperature of 26.60 ± 0.72 °C compared to 26.20 ± 1.20 °C in the SML. For salinity, significant differences were detected between SML and ULW during Cycles 1 and 2, with ULW exhibiting higher salinity levels in Cycle 1 (SML: 38.90 ± 0.32 g kg⁻¹; ULW: 39.10 ± 0.32 g kg⁻¹) and in Cycle 2 (SML: 39.2 ± 0.27 g kg⁻¹; ULW: 39.4 ± 0.29 g kg⁻¹). Notably, the salinity data showed greater variability, especially in Cycle 2, in the SML, whereas the ULW remained relatively constant. The pH_T25 data collected over the three cycles displayed considerable variability, with fluctuations observed throughout the day. The SML and ULW data showed significant differences during the first two cycles, with slightly higher values in the SML for Cycle 1 (SML: 8.030 ± 0.020; ULW: 8.020 ± 0.033) and for Cycle 2 (SML: 8.020 ± 0.020; ULW: 8.010 ± 0.024). No significant differences were observed in Cycle 3 (SML: 8.020 ± 0.032; ULW: 8.020 ± 0.032).

Figure 3 Box-and-whisker plots for the different time steps at which measurements were obtained at both the SML and ULW for temperature, salinity, and pH_T25. Each boxplot represents all measurements obtained during a single 30–45 min deployment (data collected at 30 s intervals), with four daytime and two nighttime deployments per cycle. The horizontal line within each box denotes the mean of the data, whereas the vertical line associated with each box represents the 25th and 75th percentiles (Q1 and Q3) of the data.

3.3 Diel variability and depth-related anomalies

To assess diel variability in the marine carbon system, we statistically compared the thermohaline and key carbon cycle variables, pH_T25 and pCO₂, between day and night in the SML and ULW (Fig. 4, Table 2). Additionally, we calculated ΔSML-ULW to evaluate the magnitude of vertical anomalies during diurnal and nocturnal conditions. The analysis of thermohaline variables indicated significant differences between diurnal and nocturnal data at the two depths across the three cycles (Table S2). However, the salinity data recorded from ULW during Cycle 3 did not exhibit significant differences. The observed anomalies between the SML and ULW varied across the three cycles for temperature. In Cycle 1, the SML experienced positive diurnal and negative nocturnal temperature anomalies on average (0.19 ± 0.14; −0.10 ± 0.14 °C). During Cycle 2, negative diurnal SML anomalies and positive nocturnal anomalies were observed (−0.13 ± 0.12; 0.28 ± 0.18 °C). In Cycle 3, diurnal and nocturnal SML negative anomalies were detected at −0.08 ± 0.24 and −1.34 ± 0.77 °C, respectively. Likewise, it was noted that salinity anomalies in SML were negative in Cycles 1 and 2, both for diurnal (-0.23 ± 0.17; −0.14 ± 0.19 g kg⁻¹) and nocturnal data (−0.29 ± 0.30; −0.15 ± 0.29 g kg⁻¹). These results suggest that external factors may influence thermohaline variables, affecting the pronounced temporal distribution of the diel cycle.

Figure 4 Frequency distribution of the anomaly of temperature, salinity, and pH_T25 between SML and ULW observed during the day (orange) and at night (blue) across the three cycles studied.

When comparing the diurnal and nocturnal data for the variables associated with the marine carbon system at each depth, we found significant differences in pH_T25, except for the ULW data during Cycle 3 (Table S2). However, no significant differences were observed in the pCO₂ data. The pH_T25 deltas (Table 2) revealed that during Cycle 1, lower pH_T25 values were recorded in the SML for diurnal data and higher for nocturnal data (−0.006 ± 0.018; 0.038 ± 0.029). During Cycle 2, the SML showed higher diurnal (0.011 ± 0.012) and lower nocturnal (−0.002 ± 0.012) pH_T25 values. Meanwhile, for Cycle 3, the SML indicated higher pH_T25 values for the diurnal data (0.013 ± 0.012) and lower values for the nocturnal data (0.026 ± 0.009) (Table 2). For the pCO₂ data during Cycles 1 and 3, lower pCO₂ values were recorded in the SML during both diurnal (−11 ± 56 and −22 ± 22 µatm, respectively) and nocturnal periods (−4 and −8 µatm, respectively) (Table 2). In Cycle 2, lower diurnal and higher nocturnal values were observed in the SML, with values of −12 ± 19 and 7 µatm, respectively. The large variability between the nighttime and daytime data distributions at both depths reflects the influence of several environmental drivers acting on short timescales.

Table 2 Mean anomalies of temperature, salinity, pH_T25, and pCO₂ between the two depths studied (Δ(SML-ULW)).

3.4 High-resolution variability of physicochemical parameters across diel cycles

To assess the variability of biogeochemical processes during the diel cycle, we present time series of temperature, salinity, and pH_T25 (Fig. 5). In the time series, large fluctuations were primarily observed in the SML, particularly during the day. This increased variability is consistent with the patterns described in the previous section, although the standard deviations between SML and ULW are similar overall (Table S2). The changes during the day occurred rapidly, with increases and decreases spanning 3–5 min intervals for the three parameters. More specifically, during sampling at 18:00 UTC of Cycle 2, we recorded variations over brief intervals, specifically showing changes of approximately 0.28 °C in temperature, 0.30 g kg⁻¹ in salinity, and 0.016 units in pH_T25. During the night, fluctuations were less frequent and of smaller magnitude, except for a sudden change at the end of the sampling conducted at 02:00 UTC in Cycle 2. At that point, there was a slight increase in pH_T25 at both SML and ULW, along with a decrease in salinity at both depths. This observation suggests that daytime surface heating, evaporation, and production processes likely result in changes in temperature, salinity, and pH_T25, which are less pronounced at night.

Figure 5 Time series data of temperature, salinity, and pH_T25 collected during Cycle 2: Diurnal (18:00 UTC, 12 August) and Nocturnal (02:00 UTC, 13 August) measurements.

To assess the variability observed in the SML, the sigma-t time series for Cycle 2 at 18:00 and 02:00 UTC was plotted (Fig. 6), and evaporation rates within the SML were calculated for Cycles 1 and 2 (Table 3), as salinity data required for calculations in Cycle 3 were not available. The sigma-t time series throughout the day exhibits greater variability in the SML compared to nocturnal data, when density fluctuations decrease. However, this trend of variability was not observed in the ULW. The observed patterns of temperature, salinity, pH_T25, and sigma-t align with the calculated evaporation rates. In Cycles 1 and 2, the evaporation rates peaked at 14:00 UTC (0.043 and 0.074 mm h⁻¹, respectively). They remained high during the late afternoon at 18:00 UTC (0.042 and 0.041 mm h⁻¹, respectively), coinciding with the periods of highest solar radiation and wind speed (Fig. 2). In contrast, the evaporation rate was close to zero at night. This behaviour highlights the influence of meteorological forcing on the SML during the day, underscoring the connection between evaporation and the observed variability.

Figure 6 Time series during Cycle 2 for diurnal (18:00 UTC on 12 August) and nocturnal (02:00 UTC on 13 August) sigma-t data. Sigma-t is defined as the density at a given temperature and salinity minus 1000 kg m⁻³.

To study the gas exchange between the atmosphere and the ocean, we calculated the CO₂ fluxes and k values using a wind-based parameterisation (Wanninkhof, 2014) for all three cycles and during both day and night (Fig. 7). In Cycle 1, despite the limited data, a flux of 3.64 ± 1.15 mmol cm⁻² h⁻¹ was detected during the day, while the flux was close to zero at night. A similar pattern appeared in Cycles 2 and 3, where the fluxes peaked at approximately 14:00 UTC and declined to near zero thereafter. However, the mean flux was higher in Cycle 2 (2.04 ± 3.18 mmol cm⁻² h⁻¹ ) than in Cycle 3 (0.21 ± 0.28 mmol cm⁻² h⁻¹), consistent with stronger winds. The average wind speeds during the day were 1.5, 1.9, and 0.4 m s⁻¹ for each cycle, respectively, while at night, winds dropped to nearly 0 m s⁻¹ in the three cycles. The k values during the day were 1.12 ± 0.15, 2.22 ± 3.31, and 0.09 ± 0.11 cm h⁻¹, whereas at night, they were close to 0 cm h⁻¹. In this context, excluding nocturnal fluxes in the daily average calculations introduced local percentage errors of 33 %, 50 %, and 43 % for Cycles 1, 2, and 3, respectively. The increased daytime wind speeds enhanced CO₂ fluxes, whereas calm nighttime conditions were associated with reduced or nearly zero gas transfer velocities, consistent with the observed diel SML and ULW pCO₂ variability (Fig. S2).

Figure 7 CO₂ fluxes between the atmosphere and the ocean during the studied cycles.

Table 3 Estimated evaporation rates (mm h⁻¹) based on latent heat flux and seawater density in the SML.

4 Discussion

This study revealed high variability in the SML and ULW during the diel cycle, with significant differences (p<0.05) in temperature, salinity, and pH_T25 when comparing diurnal and nocturnal data at both depths (Table S1). These results highlight the differences between the meteorological forces that influence the physicochemical properties of seawater during the day and night. During the day, the combined forcing of solar radiation and increased wind speed enhances evaporation rates (Table 3), leading to cooling the SML (Gassen et al., 2023), and driving short-term changes in temperature, salinity and pH_T25 (Fig. 4). Although temperature strongly influences in situ carbonate chemistry (Zeebe and Wolf-Gladrow, 2001), pH_T25 variations observed here do not result from a direct effect on carbonate speciation, but from temperature driven physical processes such as evaporation, which tends to increase salinity, alkalinity and thus pH_T25, and vertical mixing, which can dilute this signal by introducing water with lower pH_T25 and higher pCO₂. These patterns reflect the stronger influence of surface forcing on the SML compared to the ULW. This means that the SML is subjected to short-term fluctuations that coincide with changes in wind speed. These fluctuations influence thermohaline features, CO₂ system parameters (Acuña et al., 2008), and the kinetics of the metabolic processes occurring in the marine environment (Nimick et al., 2011). Accordingly, in response to the solar photocycle, many marine biogeochemical processes operate on a 24 h cycle (Nimick et al., 2011), with daily variations comparable in magnitude to the annual variations associated with the amount of solar radiation reaching the ocean surface at different times of the year (Herring et al., 1990; Nimick et al., 2011).

Regarding the variability observed during the day and night, we detected differences in the diurnal temperature data, reaching up to +1.89 °C in the SML and +1.36 °C in the ULW. Across the three cycles, the mean SML and ULW temperatures (Table S1) fall within the range expected for the Middle Adriatic Surface Water mass. Previous studies have reported mean summertime temperatures in the SML of 25.1 °C (Frka et al., 2009) and 27.4 ± 2.9 °C (Milinković et al., 2022). Interestingly, the negative diurnal SML anomalies in Cycles 2 and 3 differ from the typical low-wind patterns (Wurl et al., 2019), suggesting that localised mixing or evaporative cooling may have offset surface warming. In our study, strong stratification combined with weak winds likely promoted localised evaporative cooling and convection, offsetting the expected surface warming, as described in previous studies of near-surface instabilities (Cronin and Sprintall, 2001; Soloviev and Lukas, 2013). Similarly, the observed negative salinity anomaly in the SML could also be linked to these processes. While evaporation driven by intense solar radiation tends to increase salinity in the SML (Frka et al., 2009; Wurl et al., 2019), the consistently negative salinity anomalies observed in Cycles 1 and 2 (Table 2) suggest active vertical mixing, most likely associated with convective processes. This explains the daily variability in salinity distribution between the SML and ULW, with mean values of 39.09 ± 0.33 and 39.15 ± 0.34 g kg⁻¹, respectively.

In this context, the contribution of external forcing factors, such as river discharge, precipitation, and tides, to the observed diurnal variability appears to be minimal. The Krka River has the highest outflow rates in the region, but it typically experiences lower discharge in the summer. This seasonal decrease in freshwater input contributes to the higher salinity values observed in both the SML and ULW compared to periods of stronger river discharge (Frka et al., 2009; Marcinek et al., 2020). Submarine groundwater discharge is also minor, between 0.19–0.31 m³ s⁻¹ during dry periods (Liu et al., 2019), compared to its annual mean of 52.9 m³ s⁻¹ (Bužančić et al., 2016; Marcinek et al., 2020). In addition, precipitation during the first half of August 2020 was very low, averaging 1.4 mm (Croatian Meteorological and Hydrological Service, 2020), further limiting freshwater input and variability in river flow. Similarly, tidal effects at this location are minimal, characterised by a microtidal range of 0.2 to 0.5 m (Cukrov et al., 2024a). Given the stable stratification of the estuary, such small tidal amplitudes are insufficient to generate significant tidally driven vertical mixing. Instead, the observed variability can be explained by localised near-surface mixing linked to evaporation, density instabilities, and short-term turbulence at the air–sea interface. Furthermore, the observed variability follows a diurnal (24 h) cycle rather than a semi-diurnal (12 h) cycle, reinforcing the interpretation that the patterns are primarily driven by surface warming and evaporation rather than tidal forcing (Wurl et al., 2019). This contrasts with mesotidal environments (e.g., Stolle et al., 2020), where tidal mixing complicates the detection of diurnal variability, highlighting the advantage of studying these processes in microtidal settings such as the Adriatic.

The interaction between biological processes and the physicochemical properties of seawater is complex (Álvarez et al., 2014; Cantoni et al., 2012) and has a delicate balance in the marine environment (Cantoni et al., 2012; Takahashi et al., 2002). This interaction directly influences biogeochemical processes typically regulated by production and respiration (Poulson and Sullivan, 2010). These processes significantly affect the pH of seawater through the uptake or removal of CO₂. During the day, photosynthesis lowers pCO₂ levels and increases pH by consuming CO₂, whereas at night, CO₂ from respiration accumulates, increasing pCO₂ and decreasing pH (Ragazzola et al., 2021) . However, as we can see in Cycle 3 (Table S2), the observed increase in pCO₂ (525 ± 47 µatm) and pH_T25 (8.042 ± 0.020) values during the day compared to those measured at night (465 µatm; 7.993 ± 0.026) is largely consistent with the temperature dependence of CO₂ solubility, as daytime warming naturally elevates pCO₂ in surface waters (Takahashi et al., 1993; Weiss, 1974). This pattern likely reflects both the thermal effect and biological processes that can contribute to the diurnal increase in pCO₂: enhanced respiration and CO₂ accumulation in the upper layers of the water due to limited mixing may offset the anticipated photosynthetic uptake (Takahashi et al., 2002), while photoinhibition under solar radiation (Feng et al., 2008) could reduce photosynthesis efficiency. Taken together, these temperature-driven and metabolic processes provide a plausible explanation for the unexpected daytime increase in pCO₂ despite favourable light conditions (Takahashi et al., 2002).

The complexity of the coupled thermohaline and pH dynamics in seawater is highlighted in the time-series results (Fig. 5). The observed fluctuations in the SML for temperature, salinity, and pH_T25 may be due to buoyancy fluxes. Wind, thermohaline fluctuations, precipitation, and evaporation have a significant influence on surface turbulence (Cronin and Sprintall, 2001). The SML absorbs heat from sunlight and cools through radiation and heat loss, leading to changes in temperature and salinity that disrupt buoyancy, cause convective overturning, entrain deeper water from the ULW, and eventually promote mixing (Cronin and Sprintall, 2001). Wind can also enhance this process by creating tangential stress that acts as a vertical momentum flux. Temperature and salinity changes in the SML led to stratification or convection, with mixing, depending on oceanic and atmospheric forcing (Fig. 5). This process has already been observed in the SML (Wurl et al., 2019) and was found to regulate buoyancy fluxes through evaporative salinisation, playing a crucial role in the exchange of climate-relevant gases and heat between the ocean and the atmosphere.

In this context, the existence of these buoyancy fluxes only during the day could also explain the diel difference in CO₂ exchange between the atmosphere and the ocean. As observed in this study, the CO₂ fluxes exhibited differences between daytime (1.94 ± 2.45 mmol cm⁻² h⁻¹) and nighttime (0.01 ± 0.01 mmol cm⁻² h⁻¹) conditions. This pattern is consistent with the strong dependence of the gas transfer velocity (k) on wind forcing (Wanninkhof, 2014), confirming that wind is the dominant driver of diel variability. In addition, daytime buoyancy fluxes, enhanced by wind-driven evaporation and density instabilities, may further facilitate CO₂ exchange during the day, amplifying the effect of wind. Additionally, the absence of wind at night reduced the calculated flux to nearly zero. This observation is consistent with previous suggestions that gas transfer velocity parameterisations without an intercept may underestimate conditions in very calm waters (Ribas-Ribas et al., 2019). Other processes, such as small-scale convection or vertical pCO₂ gradients in the upper water column (Liss and Merlivat, 1986; Stolle et al., 2020), may also modulate short-term variability, but their contribution in our dataset appears minor compared to wind forcing. At the same time, SML properties themselves may contribute: slightly lower temperature and salinity enhance CO₂ solubility, while reduced turbulence at the boundary layer further limits exchange.

Consistent with this, the SML pCO₂ shows a considerable diel variability with a clear day–night shift in Cycle 1 (∼ 53 µatm), but smaller in Cycles 2 and 3 despite larger short-term amplitudes. Although these processes influence pCO₂, their impact on CO₂ fluxes remains small compared to the effect of k. This explains why flux differences are primarily driven by wind forcing. These patterns indicate that changes in solubility caused by temperature and the buoyancy-induced mixing determine most of the observed variability, rather than biologically driven changes. As a result, using only daytime pCO₂ values leads to a slight overestimation of daily mean ΔpCO₂ (e.g., +13 µatm in Cycle 1), but this effect is small compared to the much larger bias caused by neglecting nighttime fluxes. Temperature effects alone can account for a substantial fraction of diel pCO₂ variability, given the well-established ∼ 4 % °C⁻¹ sensitivity of seawater pCO₂ to temperature (Takahashi et al., 2002; Zeebe and Wolf-Gladrow, 2001).The main implication of our results is that neglecting nocturnal fluxes leads to a systematic bias: the three diel cycles were overestimated by 33 %, 50 %, and 44 %, respectively. Thus, while wind remains the principal driver, accounting for variability in both k and pCO₂ is essential to avoid biased daily estimates.

5 Conclusions

This study observed a clear diel variability in the distribution of thermohaline features and variables describing the marine inorganic carbon cycle, including temperature, salinity, pH_T25, and pCO₂. The SML experiences pronounced fluctuations in these parameters throughout the day, influenced by daily changes that alter near-surface stratification and mixing, mainly through evaporation-driven buoyancy fluxes and wind-induced turbulence. In addition, higher CO₂ fluxes were observed during the day, coinciding with increased wind speeds and buoyancy fluxes that enhanced the exchange of CO₂ between the two compartments. Although diel differences in pCO₂ were detectable, our results suggest that this variability is largely consistent with temperature-dependent CO₂ solubility and buoyancy-induced mixing, while a persistent biological day–night signal is not evident in this dataset. Thus, by emphasising the study of diel cycles, it has been observed that the daily variability of biogeochemical processes is in a delicate balance, making it challenging to obtain a comprehensive global understanding of marine chemistry.

To enhance understanding of these findings, implementing a comprehensive 24 h observational period is essential for accurately capturing the short-term variability that influences near-surface biogeochemistry. In recent decades, there have been technological advances in sampling equipment designed for short temporal and spatial scales. However, such progress has not been extended to nighttime observations, which remain significantly more challenging due to the greater logistical complexity and heightened safety considerations associated with operating aboard oceanographic vessels at night. In this context, it is essential to study complete diel cycles, which are crucial for understanding the global carbon budget and its associated uncertainties. Thus, generating a network of diurnal cycle data will identify the drivers of changes in marine chemistry, allowing for the assessment of the responses of marine ecosystems in the context of climate change.