A Heat Budget of the Mar Menor Lagoon, Spain

1. Introduction

The Mar Menor, located in the southeast of Spain, is the second largest lagoon in the Mediterranean Sea with a surface area of about 136 km². It is restricted from the open sea by La Manga, a 22-km long sandy barrier crossed by three tidal inlets (Figure 1). As a result of an annual evaporation rate that is 3 to 4 times greater than rainfall, it is hypersaline (42-47 ppt) in parts [1,2] which gives the lagoon a unique environmental character such that it constitutes a refuge for some endemic species. Moreover, it contains 856 ha of fringing salt marsh that is used as a stop-over by migrating shore birds. As a result it was designated a Ramsar wetland in 1994. It is shallow (a mean depth of 4.5 m) and has been highly impacted by “massive urban growth” [3]. Some sites are contaminated with heavy metals from old mining discharges [4,5] manifested as bio-accumulation (of Zn, Cd, and Pb) in bivalves and the take up of toxins by seagrass plants such as Cymodocea nodosa. The solubility and mobility of such toxic metals increases with water temperature [6]. This is particularly true in fine-grained sediment such as those of the Mar Menor. The increase in eutrophication and anoxia due to “unsuitable agricultural methods”, that introduce excess nitrates and phosphates into the system, has had significant impacts on the biosphere of the lagoon [7,8,9,10]. Groundwater discharges and flood events introduce nutrients loads from the catchment area. These events can produce algal blooms that significantly reduce the water oxygen content and can cause mass mortalities as observed in recent years (2016, 2019 and 2021). Water temperature has been shown to be a major factor leading to enhanced eutrophication in aquatic settings, while in combination with the two nutrients has far reaching effects on biological community dynamics [11]. Bio-accumulation of toxins are also enhanced by the long residence time of waters (average of 0.79 years) in the [12,13].

A surface energy budget of a coastal lagoon (Nueva lagoon) in S.E. Spain was undertaken by [14], referred to hereafter as RR06. The formulation and structure of this budget has been followed in the Mar Menor by [15], referred to hereafter as MA11, who asserted that there was limited information on the hydrodynamics of the system. More recently, 3-D numerical modelling of the hydrodynamics of the Mar Menor has significantly enhanced the understanding of structure and dynamics of the lagoon under both tidal and wind-forced conditions [13]. Whilst temperature and salinity were imposed in a 3D hydrodynamic model (SHYFEM) open boundary [16], sources and sinks of heat within the system were not included. This paper attempts to provide a simple heat box model of the Mar Menor that includes all major local sources and sinks of heat following the structure defined in the previously quoted papers. The exchange of heat across the benthic boundary layer is rarely included in models due to a lack of information on the thermal properties of the bottom sediments (a notable exception is RR06). Part of this work is to measure the temperature gradient in the topmost 0.45 m of the seabed in order to determine the heat exchange (in W/m²) of the seabed with the waters of the Mar Menor which may provide useful input to advanced hydrodynamic models such as SHYFEM.

The Mar Menor is dominated by fine grained sediments. The distribution, character and dynamics of such sediments can be greatly impacted by changes in water temperature as can benthic exchanges with the water column that impacts turbidity and water clarity [10,17,18]. Seasonal fluctuations in sea surface temperature (SST) in the Mar Menor are of the order 21^oC (10 – 31^oC) largely due to the shallow water depths [9]. These fluctuations are higher than that in the open sea. Seasonal effects of SST on biological productivity of coastal lagoons can be large [19] as is the impact of dredging of the channels connecting Mar Menor with the open sea [13]. Furthermore, the SST of coastal lagoons has been increasing at a much faster rate that the global average [20,21,22]. SST in Venice lagoon, for example, has been increasing at a rate of up to 1.75^oC/decade since 2008 [21]. It is not unreasonable to expect that the temperature trend in Mar Menor will be greater than is the case for Venice lagoon as the flushing of warm water with the open sea results in a residence time of up to 1 year [13]. An objective of the work herein is to understand the SST trends on the Mar Menor and the factors that influence these as a basis for enhanced management purposes for sustainability, as recommended by [23].

2. Materials and Methods

The factors influencing the hydrology and hydrodynamics of the Mar Menor have been defined by [2]; their definition of inputs and losses to the system has been used as a guide herein. In particular, the gross annual budgets of surface runoff, new water exchange at the inlets, ground water discharge, and evaporation are taken from this source. Other sources of information are defined below.

2.1. Thermal Probes

Hourly measurements of the temperature were made by four self-logging thermistors (Signatrol® SL52T) installed in a thermal probe that was driven into the seabed in 0.3 m of water at low tide within the shallow subtidal region of the barrier Island at Playa de Galan (Lat: 37.685^oN, Lon: 0.737^oW, see Figure 1). The thermistors were factory calibrated to an accuracy of ± 0.1^oC. Each thermistor was capable of logging temperature data hourly for about 5 months (3600 data points). Thermistors were located at a height of 0.15 m above the seabed, and at depths (h) of 0.15, 0.30 and 0.45 m below the seabed. Each thermistor was pre-programmed with date/time, and the data downloaded to ASCII-format files using company software. The mean heat flux (q_z) was estimated from the temperature gradients using Fourier’s law:

$$q_z= -k{\displaystyle \frac{\partial T}{\partial z}} \mathrm{W}/{\mathrm{m}}^2$$

(1)

where k is thermal conductivity of saturated shelly sand (assumed as 1.5 W/m.K after [24]). Thermal diffusivity (α) is related to thermal conductivity as follows:

$$\alpha = {\displaystyle \frac{k}{\rho C_p}} {\mathrm{m}}^2/\mathrm{s}$$

(2)

where ρ is the saturated bed sediment bulk density of the local well-sorted fine sand (d₅₀ = 0.23 mm) and is assumed as 1700 kg/m³, and C_p is the specific heat capacity (assumed to be 2.4 MJ/m³.K). Estimates of the temperature gradients (both during the day and night) are determined by defining the peaks and troughs of each temperature time series (that fluctuate daily). First order polynomial best fits to the peak and trough temperatures were made to define daytime and night-time values of ${\displaystyle \frac{\mathsf{\delta }T}{\mathsf{\delta }z}}$. It is assumed that the peaks and troughs in temperature are coincident in time and that the properties of the bed are constant with depth. Thermal diffusivity may also be estimated from:

$$\alpha = {\displaystyle \frac{{∆z}^2}{2\ast lag}} {\mathrm{m}}^2/\mathrm{s}$$

(3)

where Δz = 0.45 m, and the lag is the time difference between upper and lower peak (in seconds).

2.2. The Mean Monthly Time Series Data

Mean monthly SST were downloaded from global sea surface temperature data sets (HadSST1.1) compiled by the Met Office, Hadley Centre (https://www.metoffice.gov.uk/hadobs/hadisst/; [25]). A sub-set of the western Mediterranean Sea (lat: 37-38^oN; lon: 0-1^o W, see Figure 1) was extracted using purpose-written code. SST for the Mar Menor were downloaded from MODIS Aqua (11 microns, night only, 4 km resolution) using the Giovanni tool (https://oceancolor.gsfc.nasa.gov/resources/giovanni/) available from July, 2002 to the present (lat: 0.854 – 0.729 W; lon: 37.633 – 37.819 N, see Figure 1). [26] have shown that MODIS SST data is highly correlated with measurements in the Mediterranean Sea. SST data in Mar Menor prior to 2002 were backfilled to September, 1981 using the OISST blended analysis from NOAA (https://psl.noaa.gov/mddb2/makePlot.html?variableID=156646; lat: 37.25 – 37.5 N; lon: 0 – 1 W). El Niño Southern Oscillation SST anomaly based on HadISST data (ENSO3.4, https://psl.noaa.gov/gcos_wgsp/Timeseries/Nin34/) (in units of degrees Celcius) were downloaded from January, 1870 to May, 2024.

Mean monthly air temperature, wind speed, rainfall and cloud cover measured at San Javier airport (see Figure 1) were provided by AEMET (Agencia Estatal de Meteorologia) for the period 1950 – 2024. The major ocean oscillations were downloaded as follows: the Atlantic Multidecadal Oscillation (AMO, https://psl.noaa.gov/data/timeseries/AMO/); the Pacific Decadal Oscillation (PDO, https://www.ncei.noaa.gov/access/monitoring/pdo/): The North Atlantic Oscillation (NAO, https://www.ncei.noaa.gov/access/monitoring/nao/); and the Indian Ocean Dipole (IOD, https://www.cpc.ncep.noaa.gov/products/international/ocean_monitoring/IODMI/DMI_month.htm).

Gross domestic product (GDP) for Spain was digitized from data provided by World Bank Group (https://data.worldbank.org/indicator/NY.GDP.MKTP.CD?locations=ES). Sunspot activity was provided by the World Data Center (WDC) Sunspot Index and Long-term Solar Observatory (SILSO), at the Royal Observatory of Belgium, Brussels (DOI: doi.org/10.24414/qnza-ac80).

2.3. Statistical Analyses

The HadISST1.1 global monthly time series was analysed using purpose written Matlab® codes. Data from January, 1900 to December, 2024 were examined in detail for this study. A 13-point moving average was applied to all timeseries to remove seasonality. The results were subtracted from the raw data to derive the seasonal components. A 13-term Henderson filter was applied to the raw data as a means of comparison of the filtered time series [27]. The long term trends of the de-seasoned data were derived by a least square linear best fit regression. A time-series analysis of each calendar month was also undertaken following [28] and a least-square trend analysis was applied. The level of significance of the trends (p-values) were examined in Sigmaplot® and associated level of significance were reported. SST anomalies (SST-a) as defined by [29] were used in statistical analyses which were defined as the de-seasoned value minus the long-term mean. The same procedure was applied to air temperature anomaly (AT-a) and CO₂ anomaly (CO2-a) to provide de-seasoned time-series that were comparable to the oceanic indices listed above.

Mean monthly data (SST-a “the dependent variable”, NAO, ENSO, IOD, AMO, PDO, CO2-a, sunspot number, GDP, AT-a, wind speed, and incoming solar radiation “the independent variables”) for two time periods (pre and post 1980) have been analysed by cross-correlation (to examine phase relationships), by Spearman rank correlation (non-parametric analysis less influenced by outliers), and by principal component analyses in order to understand which factors appear dominant to explain the variances of air temperature and SST. Cross-correlation of variables with SST-a and AT-a was evaluated in Matlab®, other statistical tests were undertaken in Sigmaplot®. [23] show that satellite derived SST can be 0.2^o K different from waters below, and may also be influenced by dust in the atmosphere. These effects have not been considered in this study.

2.4. A Box Model for the Mar Menor

A box model of energy inputs and outputs to the Mar Menor seawater has been constructed following the work of [15] in order to assess the short term (months to years) drivers of SST. The purpose of the box model is to quantify (and balance) the sources and sinks of energy that result in the observed SST and to compare these results (at monthly time steps over 22 years, 2002 to 2024) to the daily time steps (over 3 years) in [15]. A second purpose was to see if there were any fundamental changes in the heat balance that could be detected over the 22-year simulated period. Six main sources/sinks are recognised herein following [30]: (1) sensible heat exchanges between air, ground, and seawater, (2) latent heat exchanges with the atmosphere due to evaporation and condensation, (3) incoming short and long wave solar radiation (modified by cloud cover), (4) outgoing long-wave radiation (also modified by cloud cover), (5) seawater/ground exchanges in sensible heat, and (6) tidal mixing/advection across the inlets. The sources/sinks of 1 to 4 are described in full by [31,32] together with relevant constants and coefficients. Tidal mixing across the inlets is described by [13].

The sensible heat flux (Q_air) from air to seawater is dependent on the mean temperature gradient between air and water and the wind speed over water. It is defined as follows:

$$Q_{air}= {\rho }_{air}\ast {Cp}_{air}\ast C_h\ast {xU}_{10}\ast (T_{sea}-T_{air}) \mathrm{W}/{\mathrm{m}}^2$$

(4)

where ρ_air is air density (1.22 kg/m³), Cp_air is the specific heat of air (1006 J/kg.C), C_h is the bulk transfer coefficient (1.2 * 10^-3), U₁₀ is the wind speed (m/s) at a height of 10 m, assumed in the first instance to be better represented by the friction velocity (U_*) following the method of [33], x is an empirical coefficient (wind reduction factor) used to equate predicted (1.29 m/annum) and measured evaporation rates (1.35 m/annum, [2,15], and T_sea and T_air are the mean monthly seawater and air temperatures respectively.

SST in Mar Menor is also modified by the latent heat exchange due to evaporation (Q_lt) which, according to [32] requires information on the saturation specific humidity (e_sat) and saturation air vapour pressure (q_sat):

$$Q_{lt}= l_v\ast C_e\ast C_h\ast {xU}_{10}\ast (q_{sat}-q_{air}) \mathrm{W}/{\mathrm{m}}^2$$

(5)

where l_v is the latent heat of vaporization (2.43 x 10⁶ J/kg at 20^oC and 35 psu), C_e is the evaporation coefficient (1.5 * 10^-3), and where $q_{sat}=0.622\ast {\displaystyle \frac{e_{sat}}{P-0.378\ast e_{sat}}}$ where P is atmospheric pressure (assumed as 1013 millibars or hectopascals) and:

$$e_{sat}=6.112\ast \left({\displaystyle \frac{e^{\left(17.67\ast T_{air}\right)}}{{(T}_{air}+273.5)}}\right)$$

(6)

The effect of salinity on evaporation rate has been ignored.

Incoming top of the atmosphere total solar radiation (Q_insol) is read from a mean monthly file (https://asdc.larc.nasa.gov/project/CERES/CERES_EBAF_Edition4.1). The incoming short and longwave radiation (Q_sol) is calculated as:

$$Q_{sol}=Q_{insol}\ast \left(1- {\alpha }_s\right)\ast \left(1-0.7n_c\right) \mathrm{W}/\mathrm{m}$$

(7)

where α_s is the surface albedo of semi-arid terrains such as those surrounding the Mar Menor (≈ 0.25), 0.7n_c is impact of the cloud cover fraction [31]. For present purposes, a long term monthly average of cloud cover has been used ((https://gmao.gsfc.nasa.gov/reanalysis/merra-2/).

Outgoing longwave radiation (Q_out) is derived using the modified Stefan-Bolzmann relationship:

$$Q_{out}= \in \ast \sigma \ast {\left(T_{sea}+273.5\right)}^4\ast \left(0.39-0.05{q_{sat}}^{0.5}\right)\ast 0.19\left(1-0.01{n_c}^2\right) \mathrm{W}/{\mathrm{m}}^2$$

(8)

where ϵ is the emissivity of seawater (0.98), q_sat is the vapor pressure of water (assumed constant at 23 mbar) and σ is the Stefan-Bolzmann constant (5.67 * 10^-8 W/m²/K⁴) where (K = T_sea + 273.5^oC) is the absolute water temperature. The modification to EQ 8 has been proposed by [34]. This formulation assumes cloud cover influence as defined above, as recommended by [23].

Input of heat from rivers system has been evaluated based upon annual discharge data from [2] and river water temperature is equated with ground temperature derived from the MERRA-2 project (https://gmao.gsfc.nasa.gov/reanalysis/merra-2/). This input is defined as:

$$Q_r= {\displaystyle \frac{{\rho }_{water}\ast {Cp}_{water}\ast C_h\ast V\ast (T_r-T_{gr})}{Dt\ast A_{gr}}} \mathrm{W}/{\mathrm{m}}^2$$

(9)

where V is the mean monthly volume discharge of fresh water, T_gr is the mean monthly ground temperature surrounding the Mar Menor (assumed as a constant of 18^oC within the seabed of the Mar Menor), T_r is the mean monthly river water temperature (equated with air temperature), the time step, Dt = 2630016 seconds (i.e. one month), Cp_water is the specific heat of water (3993 J/kg.^oC) and A_gr is the assumed area of the Mar Menor (1.36 * 10⁸ m²).

The energy balance is then calculated following [31], RR06 and MA11 as:

$$Q_{bal}= Q_{sol}+ Q_{air}- Q_{gr}- Q_r- Q_{out}- Q_{lt} \mathrm{W}/{\mathrm{m}}^2$$

(10)

The temperature change in SST (dT) is evaluated at monthly increments of time (per m²); the affected water mass (M_water) is estimated using a mixing depth (h) of 4.5 m [3] and is defined per unit area as:

$$dT= {\displaystyle \frac{Q_{bal}}{M_{water}\ast {Cp}_{water}}} °\mathrm{C}$$

(11)

where unit water mass is:

$$M_{water}=1\ast {\rho }_{water}\ast h \mathrm{kg}.$$

(12)

The evaporation rate (E_v, m/month) is predicted as a function of the wind speed and relative humidity:

$$E_v= {\displaystyle \frac{C_e\ast x{U}_{10}\ast (q_{sat}-q_{air})}{{\rho }_{water}}}\ast dt$$

(13)

Annual evaporation is estimated by summing the monthly evaporation rates. U₁₀ (wind speed at a height of 10 m) has been optimized in order to equate EQ 13 to the annual evaporation in the Mar Menor. This results in a wind reduction factor (x) of 1.0. This factor was initially used in EQs 4, 5 and 13. Finally, the Bowen ratio (BR, the ratio of sensible to latent heat exchanges with the sea surface) was also estimated at monthly intervals.

3. Results

A stepwise approach to the SST in the Mar Menor was used to examine the heat budget post 1980. Firstly, the temperature profile measurements in the seabed are presented (in order to define a seawater/seabed heat exchange), followed by the Spearman statistical analyses, then the cross-correlation of factors defined above, and finally the balancing of the heat budget within a box model. The results are compared to the output from RR06 and MA11.

3.1. Sea Surface Temperature

The Western Mediterranean is the open boundary for the Mar Menor and hence is an important influence on the heat budget of the lagoon. Monthly mean SST data from the W. Mediterranean Sea, extracted from the HadSST1.1 global data set, are presented from 1900 to 2024 (inclusive) in Figure 2A. Note that there is no long clear trend in SST prior to circa 1980. Post 1980 there appears to be a systematic warming. It also appears that long-term variations in SST follow trends in the AMO (Figure 2C), while shorter term peaks and troughs follow those in the NAO (Figure 2B) as suggested by [26].

The annual maxima and minima in SST of the Western Mediterranean (Figure 3) clearly show no trend in either extreme prior to 1980. Furthermore, the extremes appear remarkably consistent in time. Summer warming since 1980 shows a significant trend of 0.34^oC/decade as well as reduced scatter; winter minima also show significant warming but at a lower rate of 0.21^oC/decade. [35,36,37] have shown similar warming in their data analysis which covers the period since 1979. [26] showed a stable increase in SST of 0.55^oC/decade for the period 2005 – 2019.

The SST anomalies post 1980 for the Western Mediterranean, the Mar Menor (MODIS), and the OISST (NOAA) together with the air temperature anomaly (AT-a) for the Mar Menor are presented in Figure 4. Note that the peaks and troughs are generally consistent between records though the Mar Menor is generally cooler that either local air temperature or SST further offshore. Clearly evident in all data sets is the relative cooling presumed due to the volcanic explosions of Mount Pinatubo (1991 – 1994) and Hunga Tonga in January, 2022 [38]. The least-squares best-fit trends in SST are presented for each month of the year in Figure 5. These trends are also listed in Table 1. Notice that all months of the year show significant warming with warming greatest in the summer and least in the winter.

3.2. Seawater to Seabed Heat Exchanges in the Mar Menor

Good quality SST and seabed temperature data were recorded hourly from 10:42, 5 May 2023 to 16:42, 14 September, 2023 in a shallow (0.3 m minimum depth) sandy region of the Mar Menor (See Figure 1). The time series was dominated by diurnal fluctuations of up to 5^oC, with a steady increase in temperature punctuated by two storm events (22 – 24, May and 26, August – 5 September) referred to herein as “Event 1” and “Event 2” (see Appendices 1 and 2). Figure 6A shows the SST during this period and the peak (daytime) temperatures from the four temperature sensors. The maximum air temperature from San Javier airport is also shown. It generally falls below the Mar Menor SST and the seabed temperatures within the shallow water setting of the temperature probe with notable exceptions taking place in the two events. Peaks at each depth were not synchronous but were delayed by an average of 4.2±2.3 hours at the deepest sensor. The diurnal variation in the vertical heat flux was large (-20 to 28 W/m², negative flux is warming downwards). This variation in vertical temperature difference is illustrated for Event 1 in Figure 6B.

The figure shows that the temperature gradient inverts on a diurnal frequency illustrating the rapid response of the seabed to changes in heating. That is, the seabed provides heat to the seawater at night. The least-squares best fit of the hourly temperature gradients for the measured duration showed a mean daytime value of -5.99±2.05 ^oC/m (heating of the seabed) and 2.59±1.51^oC/m at night (cooling of the seabed). The estimated vertical heat flux was -8.99±3.08 W/m² during the day and 3.09±2.27 W/m² at night. The “dampening” depth (D) of the temperature gradients is estimated from [39]: $D={\left({\displaystyle \frac{2k\ast 86400}{2pi}}\right)}^{0.5}$ to be 2 m into the seabed. Below this depth, the temperature is assumed to be constant. The vertical heat flux for the time series shows strong diurnal oscillations and rapid changes between the peak in maximum and minimum values. Figure 7 shows results for a sub-set of observations for clarity. It appears that estimates of thermal diffusivity are unstable during periods when the heat flux is changing most (during the day). The day-time peak temperature gradient, thermal conductivity (2.01±0.9 W/m.K), and heat flux (appear to be constant except for the two rainfall events defined earlier. The daytime temperature gradient appears to be highly sensitive to rainfall and less sensitive to air temperature. This is evident in the two events recorded herein that resulted in 26.8 and 28 mm of rainfall respectively. Under these events, the temperature gradient reduced to zero with a lag of 24 hours.

Similar (but smaller) responses were measured during 2nd, 6th and 12th June, 2023 for lower levels of rainfall but similar air temperatures. A similar finding was made by RR06. Time series of the daytime and night-time peak vertical temperature gradients have been plotted against mean air temperature (San Javier airport) and daily rainfall (Appendix 2). The response and lag of the daytime temperature gradient are clear. This contrasts with the night-time gradients that appear to be influenced by events up to 4 days prior to the rainfall events. Under other times, night-time temperature gradients were remarkably constant.

The thermal conductivity of the shelly fine sand of this study site should remain constant as should estimates of thermal diffusivity. A back calculation (using EQ 2) of the mean night-time estimate of conductivity (assuming a constant heat flux of 1.0 W/m²) was 0.52±0.45 W/m.K and for the daytime (assuming a heat flux of 2.4 W/m²) was 0.54±0.65 W/m.K. This is lower than the value for marine sediments reported by Foucher et al. (2002), and close to the value for seawater. The estimate of thermal diffusivity using EQ. 2 yielded a value of 6.28 x 10^-6 m²/s. A mean value of similar sediments is 1.7 x 10^-7 m²/s [40,41]. A lower value herein could result from a lower porosity or thermal capacity of Mar Menor sediment or from a stronger effect of ground water inflow.

The mean monthly air temperature data from San Javier airport for the period 1980 – 2024 (inclusive) have been de-seasoned in a similar fashion to the SST for the Western Mediterranean. The trends for each month of the year are presented in Figure 8 and tabulated in Table 1. All but the winter months show significant (p< 0.05) warming trends that peak during the summer months (May – August).

3.3. Step-Wise Regression Analysis

The SST and air temperature data were first used to examine the trends in SST for each month of the year (Figure 5 and Figure 8 respectively). These trends (^oC/year) are tabulated in Table 1 (^oC/decade) together with the level of confidence (p-value). The summer warming evident in Figure 3 is seen clearly with strongest trends between May and October in the HadSST1.1 data. The MODIS (Mar Menor) trends are greater than the Western Mediterranean. However, the level of significance of the trends are generally low. These trends correspond well to the air temperature trends for the region which do have high levels of significance.

The anomalies in air temperature (AT-a), SST (Mar Menor, SST-MM-a), SST (Western Mediterranean, SST-HAD-a) and NOAA OISST have been used to understand changes in SST in the Mar Menor. Anomalies were chosen to remove the local effects of seasonality. The post 1980 time series are shown in Figure 4. The peaks and troughs appear to coincide reasonably well. As the Mar Menor SST data set from MODIS started in July, 2002, SST prior to that time was back-filled by regression of the OISST data with the MODIS data post 2002 (r² = 0.88, n = 538; p < 0.001).

The results of the Spearman correlation analysis of NAO, ENSO3,4, IOD, AMO, PDO, CO2-a, sunspot number, GDP, AT-a and the SST anomaly in the Western Mediterranean (HAD) against SST-a in the Mar Menor are presented in Table 2A.

3.4. The Principal Component Analysis

The PCA analysis conducted on mean monthly data from 1980 to 2024 (inclusive) showed that the strongest correlation of SST in Mar Menor was with SST in the adjacent coastal waters with the Mar Menor leading by 1 month. As shown in the heat box model (see later), this is not due to mixing across the inlets. Air temperature showed the next strongest correlation with SST_MM but with zero lag. The strong correlation to CO₂ across all lags shows that both are rising in a similar fashion. None of the other ocean indexes showed significant correlations. Curiously, sunspot activity was negatively correlated at all lags supporting the contention of Vinõs (2022) that warming takes place during periods of low sun spot activity: GDP is also inversely related which is diagnostic of the fact that industry is lower during the summers. Air temperature shows similar patterns to the SST in Mar Menor, but also includes a significant contribution of AMO (with 1 to 10 month lags) and IOD (with a 1 month lag). Air temperature appears inversely related to PDO (with a 12 month lag). There is no clear relationship of air temperature with either ENSO3,4, NAO, sunspot activity, or GDP.

Figure 9. (A) Cross correlations of key parameters with air temperature (San Javier airport) for lags of -36 to 36 months. The strongest correlation (at zero lag) is with the SST in the Western Mediterranean. This is followed by SST in the Mar Menor, CO₂, AMO and IOD. Notice there is almost no correlation with GDP. (B) Cross correlations of key parameters with SST in Mar Menor. The strongest correlation is with SST (Western Mediterranean) followed by air temperature, CO₂, and AMO. Note the inverse correlation with sunspots and GDP.

Figure 9. (A) Cross correlations of key parameters with air temperature (San Javier airport) for lags of -36 to 36 months. The strongest correlation (at zero lag) is with the SST in the Western Mediterranean. This is followed by SST in the Mar Menor, CO₂, AMO and IOD. Notice there is almost no correlation with GDP. (B) Cross correlations of key parameters with SST in Mar Menor. The strongest correlation is with SST (Western Mediterranean) followed by air temperature, CO₂, and AMO. Note the inverse correlation with sunspots and GDP.

3.5. The Heat Box Model

An estimation of the heat budget of the Mar Menor has been evaluated monthly from January, 2002 to May, 2025 inclusive. A balance in the heat budget was sought for the sources of heat (principally incoming solar radiation), outgoing heat (longwave radiation) and the internal exchanges of sensible (air/water, rain, river, and tidal exchanges) and latent heat losses through evaporation. The resulting predicted SST was compared against measured SST derived from MODIS. A first approximation was made to balance the published annual evaporation (m/annum) using EQ 13. Solar radiation at ground level was estimated as 80% of that at the top of the atmosphere. As well, estimates were made with and without the effects of clouds. The best closure with observations was found by neglecting the MERRA 2 average monthly cloud cover as the scatter was reduced by doing so. The same was true for estimates of the outgoing radiation. The resulting estimated time series of SST is shown in Figure 10A together with the difference between observed and predicted SST (Figure 10B). A histogram plot of the differences for the entire time series is shown in Figure 11. The scatter in the differences is Gaussian and yields a standard deviation of ±2.6^oC. Figure 10D and Figure 10E show the time series of the major sources and sinks of heat: incoming solar radiation, outgoing longwave radiation, latent heat flux and the different exchanges of sensible heat. The monthly means and annual values of each is given in Table 3.

) which is far too small to be detected within the simple model presented.

The influence of tidal exchanges at the three inlets, the sensible exchange through the seabed and the other sources (rain and rivers) were small and had minimal impact on the estimate of SST.

4. Discussion

[42] have examined SST trends in the north Atlantic during the 20th century, concluding …” no one should doubt that human activity is largely responsible” and that increasing droughts are likely in the Mediterranean region driven in part by the AMO [43]. The data examined herein (Figure 2) suggest that for the Western Mediterranean, at least, the long term trends are also driven by the AMO (which has been increasing in magnitude since circa 1970). The 1980 change in SST regime found herein for the W. Mediterranean is close to the 1976 regime change evident in the Pacific Ocean [38] which corresponds to a shift in ENSO and PDO patterns and the beginning of the current global warming trend. It is tempting to link the warming of the Mar Menor since 1980 to increases in CO₂ as there is a strong correlation between the two (i.e. they are both rising). However, the lack of a long term trend in SST pre 1980 (when CO₂ was also rising) suggests that regional factors may be more significant drivers. The systematic increase in magnitude of the AMO since circa 1970 together with peaking of the Eddy cycle of solar activity as proposed by [38] may be the main drivers of trends found herein. The AMO presently appears to be at its peak and is due to decrease during the next 30 years [38]. If this is so, then we could also anticipate a reduction in SST in the Mar Menor over this time period.

The PCA analysis conducted on mean monthly data from 1980 to 2024 (inclusive) showed that the strongest correlation of SST in Mar Menor was with SST in the adjacent coastal waters with the Mar Menor leading by 1 month. As shown in the box model, this is not due to mixing across the inlets. Air temperature showed the next strongest correlation with SST_MM but with zero lag. The strong correlation to CO₂ across all lags shows that both are rising in a similar fashion. None of the other ocean indexes showed significant correlations. Curiously, sunspot activity was negatively correlated at all lags supporting the contention of [38] that warming takes place during periods of low sun spot activity: GDP is also inversely related which is diagnostic of the fact that industry activity is lower during the summer. Air temperature shows similar patterns to the SST in Mar Menor, but also includes a significant contribution of AMO (with 1 to 10 month lags) and IOD (with a 1 month lag). Air temperature appears inversely related to PDO (with a 12 month lag). There is no clear relationship of air temperature with either ENSO3,4, NAO, sunspot activity or GDP.

The fluxes of energy into (during the day) and out of (during the night) the Mar Menor seabed fall within reasonable limits. The delay in peak temperature has not been considered as a factor in the estimation of these fluxes. The sensors were initiated sequentially with an approximate 5 minute delay in each case. The sample rate of 1 hour is considered too coarse to evaluate the travel time of heat with depth.

The heat box model for the seawater of the Mar Menor has been used to balance the sources and sinks of heat. Initially, an evaporation rate of 1.3 m [2] was used to back calculate the relevant mean monthly wind speed given by MERRA 2 at a height of 10 m (U₁₀) using EQs 5 and 13 to solve for the cumulative evaporation depth (h). This choice proved to predict an evaporation rate of 0.11±0.06 m/month and a mean evaporation rate of 1.29 m/year. The predicted evaporation rate varied considerably with season from a minimum of 0.036 m/year in December to a maximum of 0.20 m/year in August. The years 2002 and 2012-2014 exhibited relatively low summertime evaporation rates.

A multi-year trend in the energy balance appears to be present in the time series (Figure 10C). This seems related to changes in the summertime peak in latent heat loss (Figure 10D). Notice that years 2012 to 2014 summertime heat losses are relatively low, while 2006 to 2010 and post 2016 summertime heat losses are relatively high. The air/seawater exchange of sensible heat (Figure 10E) also appears to be increasing with time: it is low and with low seasonal variation from 2002 to 2007, and generally increasing subsequently. The histogram of the predicted differences between model output and measured SST in Figure 11 shows that a 3-month smoothing of predictions reduced the scatter. This suggests a month-long residual lag in response to heating not represented in the monthly mean data. The addition of tidal mixing results in a slight cool bias in results. However, it could be omitted from the heat budget estimate (as done by MA11) without loss of information. The overall annual loss of 0.52 W/m² in the heat budget estimate is within the estimated error margin and so it may be assumed that the model defines a thermal balance. A mean of 0.43^oC/decade warming of the Mar Menor seawater (see Table 1) for a mean water depth of 4.5 m equates to a long-term average heat input of 0.024 W/m² which is far too small to be detected within the simple model presented herein.

Rainfall has a significant impact on the seabed temperature profile. During event 1 (23 – 26 May, 2023), the daytime temperature gradient was almost eliminated whereas the seabed relative warming at night was enhanced to a maximum of 2.4^oC. The same occurred during event 2 (26 – 28 August, 2023), but to a greater degree. During this event the seabed remained consistently warmer than the seawater to a maximum of 5.0^oC. The heat is subsequently returned to the water column in the following 24 hours. The air temperature and wind speed were not unusual during these events.

[15] ignored the effects of freshwater inflow and rainfall in their estimates of the heat budget. The data presented herein (Table 3) suggest that to do so would have no material impact on their results. The impact of heat exchange with the ground, however, can have considerable impact on the day-night vertical temperature gradients at diurnal time scales if not longer; cooling seawater by day and warming it by night. MA11 estimated the mean annual Bowen ratio to be 0.19, which is typical of Mediterranean coastal lagoons: data herein suggests it varies throughout the year from -0.26 in December to 0.33 in April (mean = 0.053±0.21). MA11 show a slightly lower latent heat loss then found herein (see Table 3), which can be explained through differences in the air vapor pressure and relative humidity. Figure 4 and Figure 10 illustrate that the period 2003 to 2006 exhibited relatively mild winters and a narrower range of seasonal SST. This is due in part to the greater than average summertime latent heat losses during this period. Much of the differences between MA11 and those presented herein are considered due to the decadal changes in heat sources and sinks much of which is considered due to changes in the impact of the AMO.

[14] have estimated that the heat loss as long-wave radiation from a small lagoon in Spain was 70% of non-advective heat losses, and evaporation accounted for 20% of total energy heat losses. The values estimated herein were 26% and 70% respectively (see Table 3) i.e. the inverse. The remarkable similarity in heat fluxes presented by RR06 to those found herein, suggests similar ratios. The systematic increase in annual heat flux due to evaporation since 2006 (evident in Table 3) may be due to the long-term warming trends in air temperature and SST.

5. Conclusions

Long-term records of SST from the Western Mediterranean indicate that there was no trend prior to circa 1980 but appeared to respond on a decadal scale to the AMO and on an annual scale to the NAO [26,44]. Post 1980 SST showed significant increases of 0.34^oC/decade (maximum) and 0.21^oC/decade (minimum). This compared with 0.39^oC/decade for the Mediterranean Sea as a whole (Cusack and Cox, 2025) There were significant (p < 0.05) increases for all months of the year. The onset of warming appears to lag the 1976 inflection in SST in the Pacific Ocean and is attributed to changes in ENSO and PDO [38].

Air temperature trends in the Mar Menor region have shown significant systematic warming since 1980 at all months of the year. Peak warming takes place in the summer months (maximum of 0.61^oC/decade in July).

A spearman rank correlation of the SST anomaly in the Mar Menor (2000 - 2024) showed strong and significant correlations with (in order) SST in the Western Mediterranean, CO₂, Spanish national GDP, and the AMO. It was inversely correlated with the PDO and sunspot activity. A similar finding resulted from a PCA analysis of variables which also showed the importance of the IOD.

Measured heat exchanges between the seabed and the water column in the Mar Menor showed strong diurnal inversions of temperature gradient resulting in a daytime heat flux of -8.99±3.08 W/m² and a night-time reversal of 3.89±2.27 W/m² resulting in a warming of the seabed with time. The estimated depth of temperature fluctuations in the seabed was 2 m. The 5 months of hourly measurements showed that rainfall had the greatest effect on seabed temperature and heat exchange with the water column. Also clear, was the fact that a sampling rate of 1/hour is insufficient to accurately predict thermal diffusivity of conductivity of the seabed sediments. A reasonable estimate of heat flux, however, was achieved using the diurnal maxima and minima at a sample rate of 1 /hour.

A heat budget box model of SST in the Mar Menor was developed as an extension of that presented by MA11. Predictions using added inputs that included tidal mixing, seabed/water exchanges and river inflow. A balance to within ±2.6^oC was derived over a prediction period from January, 2002 to December, 2024 resulting in a insignificant heat loss of -0.52 W/m². Results herein agree with MA11 that the dominant sources and sinks of heat are solar radiation and evaporation; other terms appear insignificant to the outcome. Incoming solar and outgoing latent heat losses (through evaporation) were close to those estimated by RR06 in an adjacent lagoon and based upon daily measurements made over a single year. The use of mean monthly data (herein) thus appears to be justified in terms of examining the heat budget of lagoons such as the Mar Menor.