# Observation and Projection of Marine Heatwaves in the Caribbean Sea from CMIP6 Models

David Francisco Bustos Usta, Rafael Ricardo Torres Parra, Lien

Rodríguez-López, Maibelin Castillo Alvarez, Luc Bourrel

► To cite this version:

David Francisco Bustos Usta, Rafael Ricardo Torres Parra, Lien Rodríguez-López, Maibelin Castillo Alvarez, Luc Bourrel. Observation and Projection of Marine Heatwaves in the Caribbean Sea from CMIP6 Models. Remote Sensing, 2024, 16 (13), pp.2357. 10.3390/rs16132357. insu-04831470

HAL Id: insu-04831470

https://insu.hal.science/insu-04831470v1

Submitted on 11 Dec 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

Article

# Observation and Projection of Marine Heatwaves in the Caribbean Sea from CMIP6 Models

David Francisco Bustos Usta <fnref n="1" />,* , Rafael Ricardo Torres Parra <fnref n="2" />, Lien Rodríguez-López <fnref n="3" />, Maibelin Castillo Alvarez <fnref n="1" /> and Luc Bourrel <fnref n="4" />

<figure><img src="image_3.png" /></figure>

check for
updates

**Citation:** Bustos Usta, D.F.; Torres Parra, R.R.; Rodriguez-López, L.; Castillo Alvarez, M.; Bourrel, L. Observation and Projection of Marine Heatwaves in the Caribbean Sea from CMIP6 Models. *Remote Sens.* **2024**, *16*, 2357. https://doi.org/10.3390/rs16132357

Academic Editors: Jorge Vazquez and
Eileen Maturi

Received: 17 May 2024

Revised: 14 June 2024

Accepted: 25 June 2024

Published: 27 June 2024

<figure><img src="image_4.png" /></figure>

**Copyright:** © 2024 by the authors.
Licensee MDPI, Basel, Switzerland.
This article is an open access article
distributed under the terms and
conditions of the Creative Commons
Attribution (CC BY) license (https://
creativecommons.org/licenses/by/
4.0/).

**Abstract:** In recent decades, climate change has led to ocean warming, causing more frequent extreme events such as marine heatwaves (MHWs), which have been understudied in the Caribbean Sea. This study addresses this gap using 30 years of daily sea surface temperature (SST) data, complemented by projections for the 21st century from nineteen Coupled Model Intercomparison Project Phase 6 (CMIP6) models. In the 1983–2012 period, significant trends were observed in the spatially averaged MHWs frequency (1.32 annual events per decade and node) and mean duration (1.47 ± 0.29 days per decade) but not in mean intensity. In addition, MHWs show large monthly variations in these metrics, modulated by interannual and seasonal changes. MHWs seasonality is different in the three used metrics, being more intense and frequent in warm and rainy months (intensity between 1.01 to 1.11 °C, duration 6.79 to 7.13 days) and longer lasting in late boreal winter (intensity between 0.82 to 1.00 °C, duration 7.50 to 8.31 days). The MHWs behavior from two extreme months show that these events can occur in both small and large areas in the Caribbean. Overall, models tend to underestimate the annually averaged MHWs frequency and intensity, while they overestimate duration when compared to observations. MHWs projections are more extreme under SSP585, as they are sensible to the radiational scenario. However, an increase in MHWs intensity and duration (events lasting as much as 154 days by 2100) is expected, driving a decrease in frequency (–37.39 events per decade under SSP585 by 2100). These projections imply that MHWs conditions at the beginning of the century will be nearly permanent in the Caribbean’s future. Nonetheless, caution is advised in interpreting these projections due to differences between models’ simulations and observed data. While advancements in oceanic models within CMIP6 demonstrate progress compared to previous CMIP initiatives, challenges persist in accurately simulating extreme events such as marine heatwaves.

**Keywords:** sea surface temperature extremes; marine heatwaves metrics; MHW intensity; MHW frequency; MHW duration

## 1. Introduction

For over a century, the world's oceans have been warming due to climate change effects, generating more frequent extreme phenomena [1-4]. This warming has increased the marine heatwaves' (MHWs) intensity and duration [5-7]. At a global scale, MHWs events in the 1987-2016 period increased over 50% when compared to the 1925-1954 period, showing the Caribbean Sea (CAR hereinafter) as a hot-spot region due to the significant negative impacts observed in the biosphere [8].

The MHWs term was first introduced by [9], who quantitatively defined this event as a discrete period of anomalously warm water exceeding either a fixed [7], a seasonally

<footnote n="1">¹ Facultad de Oceanografía, Universidad de Concepción, Concepción 4030000, Chile; mcastilloa@udec.cl</footnote>

<footnote n="2">² Grupo de Investigación en Geociencias GEO4, Departamento de Física y Geociencias, Universidad del Norte, km 5 vía Puerto Colombia, Barranquilla 25260, Colombia; rrtorres@uninorte.edu.co</footnote>

<footnote n="3">³ Facultad de Ingeniería, Arquitectura y Diseño, Universidad San Sebastián, Lientur 1457, Concepción 4030000, Chile; lien.rodriguez@uss.cl</footnote>

<footnote n="4">⁴ Géosciences Environnement Toulouse, UMR 5563, Université de Toulouse, CNRS-IRD-OMP-CNES, 31000 Toulouse, France; luc.bourrel@ird.fr</footnote>

<footnote>* Correspondence: davidbustos@udec.cl</footnote>
varying [5], or a cumulative threshold [10]. MHWs are caused by oceanic or atmospheric processes only or by a combination of them, including geostrophic advection, mesoscale mixing, Ekman advection, turbulent heat fluxes, radiative heat fluxes, vertical mixing, and outgoing shortwave radiation, among others [6]. Moreover, MHWs may be modulated by large-scale variability [11].

Climate models are an essential tool to predict the future occurrence of MHWs as- sociated with different radiative scenarios [12]. These global general circulation models (GCMs) are crucial for forecasting future climate trends, simulating intricate interactions among ocean, atmosphere, cryosphere, and land and yielding insights into their dynamic behavior [13,14]. Thus, they are recognized as the one of the most foremost instruments for gauging the Earth’s climate response to diverse radiative scenarios, playing a pivotal role in assessing atmospheric impacts of greenhouse gases concentrations. The Coupled Model Intercomparison Project (CMIP) provides a framework for evaluating and contrasting outcomes from these models, developed by various research groups [15,16].

Several studies have shown that MHWs have increased in their frequency and days of occurrence, values which are expected to increase throughout the 21st century due to climate change effects [6,7,12,17]. However, the accurate reproduction of MHWs by climate models is challenging due to limitations in their spatial resolution, which cause an underestimation of sea surface temperature (SST) anomalies in eddy-rich regions [18] such as the CAR [19]. As a consequence of this limitation, global and regional MHWs assessment using models from the Coupled Model Intercomparison Project Phase 5 (CMIP5) tends to simulate less frequent, weaker, and longer-lasting events when compared to observa- tions [20–22]. Results from the latest CMIP6 models show similar biases to CMIP5 models, especially related to the models’ spatial resolution limitations [12]. However, CMIP6 mod- els have several improvements when compared to previous models, such as their better spatial resolution, the improved use of parameterization schemes to reproduce physical processes [16,23], and an overall better performance reproducing extreme phenomena [24].

Since the publication of CMIP6 results, various studies have assessed MHWs projec- tions at global and regional scales, showing an increase in their frequency and duration. These results seem to be related to climate change effects, as they are sensitive to the 21st shared socioeconomic pathway (SSP) used [12,25–27]. SSPs are described as stories that happen in the future and are based on five narratives describing alternative socio- economic developments, including sustainable development (SSP1), middle of the road (SSP2), regional rivalry (SSP3), inequality (SSP4), and fossil-fueled development (SSP5) [28]. However, MHWs studies using CMIP6 models in the CAR are scarce. Studying these events in the CAR is vital since MHWs can generate irreversible damage to marine life, causing severe alterations to marine biodiversity and to the composition and functioning of crucial ecosystems in this region [29,30]. For this reason, the few available MHWs studies in the CAR focus on their effects on coral bleaching and loss [31,32] and ecosystems impacts [33,34]. Conversely, the MHWs characteristics in the CAR such as frequency, in- tensity, and days of duration and their future projections with CMIP models have been poorly investigated.

The CAR is a semi-enclosed tropical body of water located east of Central America, north of South America and west of Antilles, which separates it from the Atlantic Ocean [35]. It is considered one of the largest seas in the world, with a surface area of approximately 2.75 million km². The depth of the CAR varies from shallow waters near the coast to depths of more than 7000 m in some areas (see Figure 1). This sea also has various important ecosystems, such as coral reefs, mangroves, coastal forests, and submarine plains [36]. The CAR tropical climate is characterized by warm temperatures throughout the year, with two main seasons modulated by the Intertropical Convergence Zone (ICTZ) meridional position. A dry season from December to March, which is relatively colder as it coincides with boreal winter, is characterized by strong winds and low precipitation. A rainy season from August to November has weak winds and is relatively warmer, as it coincides with boreal summer. A ~2 °C difference is observed in temperature at the mixed layer depth
between these two seasons [37]. During the rainy-warm season, the region is prone to hurricanes that significantly affect the region's ecology, economy, and culture.

<figure><img src="image_5.png" /><figcaption>**Figure 1.** Caribbean Sea bathymetry from [38]. Red polygon indicates Caribbean Sea limits used to assess marine heat waves (MHW). Boxes 1, 2, and 3 are used to represent the Venezuelan, Cayman, and Colombian Basins, respectively.</figcaption></figure>

In this study, we pursued three main objectives outlined as follows:

* Characterize MHWs in the CAR using different metrics (mean intensity, duration, and frequency) calculated from daily SST time series from 1983-2012 period;

* Evaluate the reliability of 19 CMIP6 models to simulate MHWs in the CAR during the period of SST observations;

* Assess MWH projections from 2041 to 2100 using CMIP6 models under two SSPs scenarios (SSP245 and SSP585), providing insights into future MHW variations in the CAR.

The main motivation of this paper is to improve the understanding MHWs in the CAR from 30 years of observations. However, to make this information valuable to stakeholders and enhance sustainable development and management measures, we need to understand how MHWs are projected to behave by the end of the century as a consequence of global warming. If this is achieved, actions can be taken to minimize negative consequences to the biosphere and coastal communities related to MHWs events. GCMs models are the best available tool to fulfill this goal. For this reason, we first evaluate the model's accuracy against observations and then present the models' MHWs projections. Therefore, results from this study contribute to fill the gap between the assessment of observed MHWs in the CAR and MHWs projections for the rest of the century under different global warming scenarios, which are an important asset to continue mitigation and adaptation efforts.

## 2. Materials and Methods

### 2.1. MHWs Definition

Fixed [7,20] or mobile thresholds [39] through the year are used to detect MHWs. Their main difference is that mobile thresholds can evolve in time, accounting for global warming, but also are used to remove the seasonality in the MHWs estimation, using a fixed reference period [7]. Therefore, the method selection will depend on the scientific question being assessed. In this study, we used mobile thresholds to remove seasonality but using a fixed period of time to better assess the MHWs forcing in the CAR, regardless of the season.
To identify MHWs, we used the method proposed by [5] to facilitate comparisons with other studies. MHWs occurs when daily SST exceeds a certain threshold for a minimum of five consecutive days. If consecutive events are interrupted for up to two days, they are considered a single event. At each grid cell, the 90th SST percentile threshold was determined using data available from an 11-day window centered on each calendar day, which was finally smoothed with a 31-day moving average. This methodology is the same as the one proposed by [22]. It was shown that modifying the thresholds percentiles does not significantly affect the conclusions regarding MHWs future projections [40].

MHWs discrete events are defined in each node by its duration and intensity. The event duration (D) indicates the number of days the SST exceeds the daily varying threshold value (Equation (1)), while the event intensity ($i_{mean}$) corresponds to the average of the daily SST anomalies through the event (Equation (4)).

$$D = t_e - t_s \qquad (1)$$

where $t_s, t_e$ represent the dates on which an MHW starts (s) and ends (e), defined in Equations (2) and (3):

$$t_s = \text{time } t, \text{ where } T(t) > T_{90}(j) \text{ and } T(t-1) < T_{90}(j) \qquad (2)$$

$$t_e = \text{time } t, \text{ where } t_e > t_s \text{ and } T(t) < T_{90}(j) \text{ and } T(t-1) > T_{90}(j) \qquad (3)$$

where $T_{90}(j) = P_{90}(X)$, in which $P_{90}$ is the 90th percentile, and $P_{90}(X)$, where $X = \{T(y,d)|y_s \le y \le y_e, j-5 \le d \le j+5\}$, with $T(y,d)$ being the SST for the year y and day of year d.

$$i_{mean} = T(t) - T_m(j) \qquad (4)$$

where $T(t)$ is the SST at the time $t$ for a given location, and $T_m(j)$ is the climatological mean, calculated over a reference period to which all the values are relative as defined in Equation (5):

$$T_m(j) = \sum_{y=y_s}^{y_e} \sum_{d=j-5}^{j+5} \frac{T(y, d)}{11(y_e - y_s + 1)} \qquad (5)$$

where $j$ is the day of year, and $y_s$ and $y_e$ are the start and end of the climatological base period, respectively.

In addition, MHWs were characterized over a given period (e.g., one year or month) in terms of the number of discrete events (frequency) as well as the mean duration and mean intensity, using a simple average of all events' duration and intensity in the period observed for each node. Spatially averaged MHWs metrics for the CAR were calculated from the simple average of the MHWs intensity and duration using all nodes in the region where at least one event occurs in the given period. Nodes without events were excluded to avoid underestimation. In the case of the spatially averaged MHWs frequency, we report the sum of events from all available nodes in the study area and period assessed [21].

Time series trends were calculated using the Theil–Sen Slope estimator [41] to avoid biases and variance in the estimation using Equations (6) and (7), where $x_j$ and $x_k$ are the data values at times $j$ and $k$ with $j > k$, $N$ is the number of samples, $i$ represents each one of the samples, and $Q_{med}$ is the final estimator. The confidence interval selected is 99% and the p-value to determine the significance was calculated using the Mann–Kendall test [42,43] defined in Equation (8), where $E[S]$ and $VAR(S)$ represent the expected value and the variance (Equations (9) and (10)), respectively, of the statistic S (Equation (11)). The advantages of these non-parametric tests include their robustness to outliers and are recommended by the World Meteorological Organization (WMO) as standard tools for trend analysis [44].

$$Q_i = \frac{x_j - x_k}{j - k} \quad \text{for } i = 1, \dots, N \qquad (6)$$
$$Q_{med} = \begin{cases} \frac{Q_{[(N+1)/2]}}{\frac{Q_{[N/2]} + Q_{[(N+2)/2]}}{2}} & \text{if N is odd} \\ \frac{Q_{[N/2]} + Q_{[(N+2)/2]}}{2} & \text{if N is even} \end{cases} \qquad (7)$$

$$Z_{MK} = \begin{cases} \frac{E[S]-1}{\sqrt{VAR(S)}} & E[S] > 0 \\ 0 & E[S] = 0 \\ \frac{E[S]+1}{\sqrt{VAR(S)}} & E[S] < 0 \end{cases} \qquad (8)$$

$$E[S] = \sum_{j=1}^{n-1} \sum_{k=j+1}^{n} \operatorname{sgn}(x_j - x_k) \qquad (9)$$

$$VAR(S) = \frac{1}{18} \left( n(n-1)(2n+5) - \sum_{k=1}^{p} q_k (q_k - 1)(2q_k + 5) \right) \qquad (10)$$

$$S = \operatorname{sgn}(x_j - x_k) = \begin{cases} 1, & x_j - x_k > 0 \\ 0, & x_j - x_k = 0 \\ -1, & x_j - x_k < 0 \end{cases} \qquad (11)$$

## 2.2. Models and Observational Datasets

We used the National Oceanic and Atmospheric Administration Optimum Interpolation Sea Surface Temperature V2 (NOAA OISSTV2) high-resolution dataset [45], to assess the current climate state using a 30-year time series (1983–2012) based on SST observations. This daily dataset, with a 0.25° spatial resolution, has been used in different studies to analyze the MHWs behavior around the planet [46–48].

Nineteen CMIP6 models from the OMIP (Ocean Model Intercomparison Project) experiment (Table 1) were used [16], considering only data from the r1i1p1f1 run. These models have the same realization, initialization, physics, and forcing. The data were accessed from the Program for Climate Model Diagnosis and Intercomparison (PCMDI) at https://pcmdi.llnl.gov/CMIP6/ (accessed on 25 October 2023). The model’s daily SSTs were downloaded from the Earth System Grid Federation (https://esgf-node.llnl.gov/projects/cmip6/, last accessed: 10 April 2023).

<table><caption>**Table 1.** List of CMIP6 models used in the study.</caption><thead><tr><th>No.</th><th>Model Acronym</th><th>Country</th><th>Horizontal Resolution (Lon × Lat in Degrees)</th><th>References</th></tr></thead><tbody><tr><td>1</td><td>ACCESS-CM2</td><td>Australia</td><td>1.0 × 1.0</td><td>[48]</td></tr><tr><td>2</td><td>ACCESS-ESM1-5</td><td>Australia</td><td>1.0 × 1.0</td><td>[49]</td></tr><tr><td>3</td><td>BCC-CSM2-MR</td><td>China</td><td>1.0 × 1.0</td><td>[50]</td></tr><tr><td>4</td><td>CanESM5</td><td>Canada</td><td>1.0 × 0.5</td><td>[51]</td></tr><tr><td>5</td><td>CESM2</td><td>USA</td><td>1.1 × 0.6</td><td>[52]</td></tr><tr><td>6</td><td>CMCC-CM2-SR5</td><td>Italy</td><td>0.9 × 0.3</td><td>[53]</td></tr><tr><td>7</td><td>CMCC-ESM2</td><td>Italy</td><td>1.0 × 0.8</td><td>[54]</td></tr><tr><td>8</td><td>EC-Earth3</td><td>UK</td><td>1.0 × 1.0</td><td>[55]</td></tr><tr><td>9</td><td>EC-Earth3-CC</td><td>UK</td><td>1.0 × 1.0</td><td>[55]</td></tr><tr><td>10</td><td>EC-Earth3-Veg</td><td>UK</td><td>1.0 × 0.8</td><td>[55]</td></tr><tr><td>11</td><td>EC-Earth3-Veg-LR</td><td>UK</td><td>0.9 × 0.9</td><td>[55]</td></tr><tr><td>12</td><td>GFDL-CM4</td><td>USA</td><td>1.0 × 1.0</td><td>[56]</td></tr><tr><td>13</td><td>GFDL-ESM4</td><td>USA</td><td>0.5 × 0.5</td><td>[57]</td></tr><tr><td>14</td><td>KIOST-ESM</td><td>South Korea</td><td>0.5 × 0.8</td><td>[44]</td></tr><tr><td>15</td><td>MIROC6</td><td>Japan</td><td>1.0 × 0.7</td><td>[58]</td></tr><tr><td>16</td><td>MPI-ESM1-2-HR</td><td>Germany</td><td>0.4 × 0.5</td><td>[59]</td></tr><tr><td>17</td><td>MPI-ESM1-2-LR</td><td>Germany</td><td>1.4 × 0.8</td><td>[60]</td></tr><tr><td>18</td><td>NorESM2-LM</td><td>Norway</td><td>1.0 × 0.5</td><td>[61]</td></tr><tr><td>19</td><td>NorESM2-MM</td><td>Norway</td><td>1.0 × 0.8</td><td>[62]</td></tr></tbody></table>
A total of 195 grid cells or nodes inside the red polygon shown in Figure 1 were used to assess MHWs from the OISST dataset. Only nodes in the CMIP6 models with >50% of oceanic area were included [49]. This information is available in the “land_area_fraction” variable provided in all CMIP6 models [50]. OISST nodes data were derived from the CMIP6 models’ reference grid.

The modeled MHWs analysis was performed using the historical run and future projections in the 2041–2100 period, which was divided in two 30-year periods to represent the mid-century (2041–2070) and end-century (2071–2100) scenarios. Projections were assessed from the two future socio-economic pathways (SSPs) SSP2-4.5 and SSP5-8.5, indicating the representative concentration pathways (RCP) 4.5 and 8.5, respectively, for the 2015–2100 period [51]. The 1983–2012 period from the historical run was used to assess the models’ performance in comparison to the OISSTV2 dataset and to compute the historical baseline climatology and the daily SST thresholds used to define MHWs in the modeled data. The MHWs assessment in 30-year periods helps to reduce the internal models’ variability and short-term fluctuations [5]. In addition, the CMIP6 models’ ability to reproduce observed SST seasonality was assessed using regional means in three areas with the same number of nodes (Figure 1). The seasonality assessment was included due to its importance in defining MHWs metrics in CAR. A regression with mean, trend, annual, and semi-annual cycles (e.g., [52]), was fitted to the spatially averaged SST time series from the OISST monthly observations and individual models in each area. These areas were selected because they have different modes of variability related to different dynamic processes (e.g., [37]).

Both the NOAA OISSTV2 and the models’ SST daily datasets were re-gridded to a common 1.0° × 1.0° grid using bilinear interpolation, avoiding extrapolation. This allows direct comparisons between both datasets. This interpolation method has been widely used in rectangular grids because it is easy to implement and reduces computational costs [53]. The 19 models in the common grid were used to compute a multi-model mean (MMM). In addition, the five models with the best performance were selected to compute their mean (5MMM). Averaging results from different individual models is a common practice that helps to reduce the inherent variability and uncertainty in predictions [16].

Models’ degree of correspondence to the NOAA OISST V2 SST dataset was evaluated with a Taylor diagram. This method uses three statistics, namely the correlation coefficient (CORR; Equation (12)), normalized standard deviation (STD; Equation (13)), and root mean square error (RMSE; Equation (14)) between models and the reference dataset. Taylor diagrams [54] provide a graphical summary of how closely a given simulated variable matches a reference one, and they are a standard method to assess models’ performance reproducing oceanographic variables such as SST [55,56] and precipitation [57], and they are also an accepted method to assess performance metrics in climate models [58,59].

$$CORR = \frac{\sum_{i=1}^{n} (y_i - \bar{y})(\hat{y}_i - \bar{y})}{\sqrt{\sum_{i=1}^{n} (y_i - \bar{y})^2 \sum_{i=1}^{n} (\hat{y}_i - \bar{y})^2}} \quad (12)$$

$$STD = \sqrt{\frac{1}{n-1} \sum_{i=1}^{n} (y_i - \bar{y})^2} \quad (13)$$

$$RMSE = \sqrt{\frac{1}{n} \sum_{i=1}^{n} (y_i - \hat{y}_i)^2} \quad (14)$$

where $y_i$ and $\hat{y}_i$ represent the observed and model values, respectively; n is the number of observations and $\bar{y}$ the mean.

In this study, seasonal assessment was performed using the average between December to March (DJF) and September to October (SON) for the dry and rainy seasons, respectively [60].
## 3. Results

### 3.1. Evaluation of CMIP6 Models Performance

CMIP6 models' performance reproducing daily SST in the CAR was evaluated against the NOAA OISST V2 using a Taylor diagram (Figure 2, Table 2). Most models showed good correspondence since statistics ranges were 0.50–0.80 for correlation, 0.61–1.27 °C for RMSE, and 0.81–1.79 °C for normalized standard deviation. The best performance was achieved by the MMM (and CESM2 model), with a correlation of 0.77 (0.80), RMSE of 0.64 (0.61) °C, and normalized standard deviation of 0.95 (0.86) °C. The CESM2 model showed a good performance reproducing MHWs in other regions, such as the South Atlantic [61], South and North Pacific, and equatorial band [62]. Conversely, the MIROC6 model has the poorer performance, with higher normalized standard deviation (1.79) and RMSE (1.27).

<figure><img src="image_6.png" /><figcaption>**Figure 2.** Taylor diagram comparing daily SST time series from the NOAA OISSTV2 dataset (red circle) and 19 GCMs (Table 1) in the 1983–2012 period. Models are differentiated by symbol and color. The 5 models' mean (5MMM—red circle) and 19 models' mean (MMM—orange square) are highlighted. Blue dotted lines represent the 0.70 and 0.90 correlation values, and gray concentric circles represent the root mean square error (RMSE). The normalized standard deviation is shown with arcs of circles connecting values in the horizontal and vertical axes. Statistics' values are presented in Table 2.</figcaption></figure>
<table><caption>**Table 2.** Normalized standard deviation (STD in °C), correlation (CORR), and root mean square error (RMSE in °C) estimated for CMIP6 models against NOAA OISST V2 dataset in the 1983–2012 period. * indicates models selected to compute the five models’ mean (5MMM).</caption><thead><tr><th>Model</th><th>STD</th><th>CORR</th><th>RMSE</th></tr></thead><tbody><tr><td>ACCESS-CM2</td><td>1.10</td><td>0.50</td><td>1.07</td></tr><tr><td>* ACCESS-ESM1-5</td><td>1.05</td><td>0.74</td><td>0.74</td></tr><tr><td>BCC-CSM2-MR</td><td>1.38</td><td>0.78</td><td>0.84</td></tr><tr><td>CanESM5</td><td>1.60</td><td>0.74</td><td>1.12</td></tr><tr><td>* CESM2</td><td>0.86</td><td>0.80</td><td>0.61</td></tr><tr><td>CMCC-CM2-SR5</td><td>1.59</td><td>0.76</td><td>1.07</td></tr><tr><td>CMCC-ESM2</td><td>1.27</td><td>0.77</td><td>0.81</td></tr><tr><td>EC-Earth3</td><td>1.46</td><td>0.76</td><td>0.93</td></tr><tr><td>EC-Earth3-CC</td><td>1.39</td><td>0.76</td><td>0.85</td></tr><tr><td>EC-Earth3-Veg</td><td>1.47</td><td>0.75</td><td>0.91</td></tr><tr><td>EC-Earth3-Veg-LR</td><td>1.46</td><td>0.76</td><td>0.90</td></tr><tr><td>GFDL-CM4</td><td>1.11</td><td>0.73</td><td>0.92</td></tr><tr><td>* GFDL-ESM4</td><td>0.98</td><td>0.72</td><td>0.76</td></tr><tr><td>KIOST-ESM</td><td>1.42</td><td>0.70</td><td>1.02</td></tr><tr><td>MIROC6</td><td>1.79</td><td>0.73</td><td>1.27</td></tr><tr><td>MPI-ESM1-2-HR</td><td>1.21</td><td>0.72</td><td>0.83</td></tr><tr><td>MPI-ESM1-2-LR</td><td>1.52</td><td>0.74</td><td>1.06</td></tr><tr><td>* NorESM2-LM</td><td>0.81</td><td>0.68</td><td>0.76</td></tr><tr><td>* NorESM2-MM</td><td>0.84</td><td>0.69</td><td>0.74</td></tr><tr><td>MMM</td><td>0.95</td><td>0.77</td><td>0.64</td></tr><tr><td>5-MMM</td><td>0.74</td><td>0.81</td><td>0.59</td></tr></tbody></table>

Moreover, we used the five models with the best performance (CESM2, NorESM2-LM, NorESM2-MM, GFDL-ESM4, and ACCESS-ESM1-5) to compute their SST mean (5MMM). The 5MMM has the highest correlation and lowest RMSE (Table 2) out of all models, improving in two of the three statistics when compared to the MMM. Therefore, hereinafter, we used the 5MMM to assess MHWs from the CMIP6 models. In addition, we assessed individual results from the CESM2 since it demonstrated good performance in nearby regions.

Furthermore, upon a reviewer's recommendation, we assessed the ability of the CMIP6 models to reproduce the SST seasonality in the 1983–2012 period in three areas shown in Figure 1. Time series from the five models with best performance were assessed (Table 2). The first area shows seasonality related to Venezuela Basin, the second area the Cayman Basin, and the third area the southern Colombian Basin. The findings reveal that both the annual and semi-annual cycles play significant roles, explaining 78% to 91% of the daily SST variance in the three distinct areas (Table S1, Supplementary Material). Notably, the annual cycle consistently outweighs the semi-annual cycle across these regions, peaking in September. There are notable differences among regions; for instance, Cayman Basin shows the highest annual cycle but a lower semi-annual cycle, while Colombia Basin shows the opposite pattern. Models generally accurately capture these cycles, although slight differences exist in their amplitude, phase, and explained variance. Overall, the models effectively simulated the seasonal SST cycle within these studied areas.

## 3.2. MHWs from Observations in the 1983–2012 Period and Comparison to Models

### 3.2.1. Annual Variability

The CAR NOAA OISSTV2 annual time series reveals positive and significant trends for the MHWs spatially averaged duration and frequency (sum of events per year from all CAR nodes) in the 1983–2012 period (Figure 3). The frequency of occurrence increases at a rate of 259.20 ± 72.93 annual events per decade (1.32 annual events per decade for each of the 195 nodes) and the mean duration at 1.47 ± 0.29 days per decade. MHWs frequency appears more influenced by interannual variability, for example, in 2005, when 1525 MHWs events occurred, considerably above the mean of 362.77 events per year
(Table 3). However, a not statistically significant increase in the mean intensity of these events was observed, which was unexpected, as significant and positive SST trends in the Caribbean were reported for the 1960–2005 and 2005–2050 periods due to climate warming [60]. As suggested by a reviewer, we investigated the possible reasons behind the lack of a significant trend in the MHWs mean intensity. We found small differences between the probability density functions computed for the first and second 15-year periods of our time series using only SST values that were part of an MHW event. These results indicate small temporal differences in extreme SST data.

<figure><img src="image_7.png" /><figcaption>**Figure 3.** Caribbean Sea annual time series of MHWs properties: (A) spatially averaged mean intensity, (B) frequency (sum of events per year from all Caribbean nodes), and (C) spatially averaged mean duration. Time series from CESM2 model (dark red), 5MMM (gray), and for NOAA OISST V2 (black). Trends estimated over the 1983–2012 period are shown in each panel for each time series in units/decade. * indicates significant trends at the 99% level.</figcaption></figure>
<table><caption>**Table 3.** MHWs properties for the 1983–2012 period in the Caribbean Sea from the OISSTV2, 5MMM, and CESM2 annual time series. Intensity and duration from spatially averaged annual values. Frequency from the sum of events per year from all Caribbean nodes. Metrics include the 30-year mean, standard deviation (Std), maximum (MaV) and minimum (MiV) values, and their year of occurrence.</caption><thead><tr><th rowspan="2"></th><th colspan="3">Intensity (°C)</th><th colspan="3">Frequency (Total Events)</th><th colspan="3">Duration (Days)</th></tr><tr><th>OISST</th><th>5MMM</th><th>CESM2</th><th>OISST</th><th>5MMM</th><th>CESM2</th><th>OISST</th><th>5MMM</th><th>CESM2</th></tr></thead><tbody><tr><td>Mean</td><td>0.94</td><td>0.41</td><td>0.59</td><td>362.77</td><td>244.01</td><td>208.73</td><td>7.96</td><td>15.97</td><td>14.71</td></tr><tr><td>Std</td><td>0.06</td><td>0.02</td><td>0.06</td><td>409.29</td><td>291.19</td><td>259.84</td><td>1.88</td><td>7.97</td><td>6.48</td></tr><tr><td>MaV</td><td>1.08</td><td>0.45</td><td>0.80</td><td>1525</td><td>929</td><td>888</td><td>14.25</td><td>38.50</td><td>31.35</td></tr><tr><td>Year MaV</td><td>2005</td><td>1986</td><td>2007</td><td>2005</td><td>2004</td><td>1988</td><td>2010</td><td>2005</td><td>2008</td></tr><tr><td>MiV</td><td>0.83</td><td>0.37</td><td>0.51</td><td>7.00</td><td>0.00</td><td>0.00</td><td>5.55</td><td>5.00</td><td>5.95</td></tr><tr><td>Year MiV</td><td>2009</td><td>2006</td><td>2009</td><td>1991</td><td>1985</td><td>1986</td><td>1989</td><td>1999</td><td>1998</td></tr></tbody></table>

We assessed the CMIP6 model's (historical run) ability to reproduce the CAR spatially averaged MHWs properties observed in the same period. Trends from the 5MMM ensemble and CESM2 model indicate an increase in the MHWs frequency and duration of 53.70 ± 51.84 to 203.63 ± 49.26 events and 0.92 ± 0.83 to 3.08 ± 1.61 days per decade, respectively. Furthermore, mean intensity trends were also insignificant (0.002 ± 0.003 and 0.005 ± 0.014 °C per decade), as found in the OISST V2 dataset. Therefore, the 5MMM and CESM intensity, duration, and frequency trends were similar to the observed values and not statistically different (except for the frequency trend between OISST and CESM2).

The models show a smaller mean in the number of annual MHWs (208.73 and 244.01 events per year) when compared to observations (362.77 events per year). Similarly, the models-averaged intensity (0.41 and 0.59 °C) is about half the observed mean (0.94 °C). In contrast, the MHWs mean duration in the models (15.97 and 14.71 days) is about twice the mean observed duration (7.96 days). The 5MMM and CESM2 time series show a higher interannual variability in the annual MHWs duration (Figure 3), as shown by the comparison between the model's range (33.5 and 25.4 days, respectively) and the OISST range (8.70 days in Table 3). The opposite is observed in the frequency, where the model's interannual variability (range between maximum and minimum yearly values) is smaller than in the observed frequency range. Results from other models used to construct the 5MMM can be observed in Table S2 from the Supplementary Material.

### 3.2.2. Interannual and Sub Annual Variability

With the spatially averaged monthly time series from OISST, we assessed the interannual variability and seasonality of the MHWs metrics in the CAR for the 1983–2012 period (Figure 4). MHWs are frequent in the CAR throughout most months; however, there are also months and longer periods (e.g., 9 months in 1991) where they are absent.

Interannual differences in the MHW metrics in the CAR were assessed using two 15-year periods: 1983–1997 and 1998–2012. Results indicate small differences in the monthly mean intensity with values of 0.97 and 1.00 °C, respectively. In contrast, the monthly mean duration and number of events in the 1998–2012 period are 8.08 days and 53.53 events per month, which is higher than in the former period, with differences of 1.68 days and 42.86 events per month, respectively (Figure 4). These results are coherent with the trends reported in Figure 3. Furthermore, 2010 was a year with many MHWs of long duration in the CAR (Figure 4), which also can be seen in the metrics' annual means (Figure 3). A similar behavior was seen in 2005, a warm year in this region with severe reported impacts in coral bleaching [10] and the hurricane season [63].
<figure><img src="image_8.png" /><figcaption>**Figure 4.** Monthly MHWs metrics in the Caribbean Sea from the NOAA OISST V2 dataset to assess inter-annual and sub annual variability. Columns represent mean intensity (°C) (A,D), frequency (events per month) (B,E), and duration (days) (C,F), respectively. First row shows the monthly means (A–C). Second row shows monthly metrics distribution in the 1983–2012 period (D–F).</figcaption></figure>

Seasonality is different in the three used metrics. Months with maximum (minimum) mean intensity are between May and September (October and April), with values ranging from 1.01 to 1.11 °C (0.82 to 1.00 °C), with the maximum (minimum) peak in June (January) (Figure 4A). Therefore, MHWs intensity is greater in the rainy-warm season, while it is smaller in the dry-cold season. July's intensity reduction compared to June and August might be associated to the activation of the Mid-Summer Draft, whose stronger winds enhances coastal upwelling and vertical mixing, cooling its surface [64]. Monthly averaged MHWs intensity can be >2 °C, as seen in February and April in 1985.

The MHWs mean frequency has higher values in two periods: between May and June, with a maximum peak of 38 events per month, and in the rainy-warm season from September to November. Conversely, lower frequency occurs in January during the dry-cold season (Figure 4B). The highest monthly averaged MHWs frequency was 328 events per month (May of 2005).

The MHWs duration and frequency can be mutually constrained, as the finding of months with a large number of events implies that the daily duration of these events is limited. The largest monthly averaged duration of an MHW (March 2010) was 20.72 days; conversely, in this month, the number of events was 171. Note that this year had a large number of MHWs events with long duration, as also seen in Figure 3. The MHWs duration seasonality shows maximum (minimum) values between December and March (May and August), with values ranging from 7.50 to 8.31 days (6.79 to 7.13 days), with the maximum peak in March (Figure 4C).

The MHWs metrics in the CAR do not show a strong seasonality, as mean monthly values ranged from 0.29 °C in intensity, 13.8 events per month in frequency, and 1.52 days in duration (Figure 4 first row). This small seasonality is related to the methodology applied in this paper to define MHWs (Section 2.1), where the threshold used to define these events varies seasonally. Still, MHWs seem to be more intense and frequent during the warm-rainy season in the CAR, while they show a larger duration at the end of boreal winter.
### 3.2.3. Spatial Variability

Additionally, a spatial assessment of the MHWs metrics was conducted during the baseline period (1983–2012), using monthly means computed from the OISST database. This is important to better understand regional oceanic and atmospheric processes that might affect sea surface temperature (SST) and trigger MHWs events in the CAR. This spatial assessment was performed using each node’s annual mean but also using the means of three months from the two dominant seasons to evaluate variations from the annual average.

The spatial average and standard deviation of MHWs intensity computed from all nodes' annual mean is 1.00 ± 0.12 °C, with nodes' values between 0.86 and 1.65 °C. The dry season (DJF) spatial average is smaller (0.84 ± 0.13 °C), with a range between 0.67 to 1.41 °C, while the rainy season (SON) spatial average is slightly greater (1.04 ± 0.10 °C), with a range between 0.88 to 1.72 °C. Therefore, spatial differences in the nodes' mean intensity in all cases is ~0.8 °C (Figure 5, first row). In addition, more intense events (> 1.32 °C) are seen in the cooler Guajira upwelling region, with larger spatial differences in the dry season (DJF). We speculate that in all cases shown, MHWs in this area would be associated with events when the wind stress and Ekman transport offshore relaxes, weakening the upwelling of relative colder waters, increasing SST enough to produce an MHW.

<figure><img src="image_9.png" /><figcaption>**Figure 5.** The MHWs spatial behavior from the NOAA OISST V2 dataset in the CAR. Columns indicate the annual, dry (DJF), and rainy (SON) seasons assessed from each node's mean in the 1983–2012 period. Rows indicate mean intensity (°C), frequency (events per month), and duration (days).</figcaption></figure>

The spatial average and standard deviation of MHWs frequency computed from all nodes' annual mean is 1.86 ± 0.23 events per month, with mean values between 1.17 and 2.33. Seasonality shows a smaller spatial average in the number of MHWs per month in the dry season (1.76 ± 0.33), while in the rainy season, a larger spatial average in the number of MHWs per month occurs (2.24 ± 0.36). In the three cases (Figure 5, second row), the spatial range among nodes is 1.16, 1.86, and 1.87 events per month, respectively, exhibiting a weak zonal pattern with a lower number of events toward the East. Thus, a smaller number of MHWs occur in the Granada and Venezuela Basins when compared to the rest of the Caribbean Sea.

Conversely, the MHWs duration (Figure 5, third row), shows a larger spatial average and standard deviation when computed from all nodes' annual mean (11.31 ± 2.11 days) and dry
season mean (10.44 ± 2.39 days) when compared to the rainy season (9.90 ± 1.62 days). The MHWs duration spatial range among nodes is largest in the annual mean (6.56 to 20.17 days) than in the two seasons assessed, which have a similar range (~11 days). A weak spatial zonal pattern can be seen, with longer-lasting events to the East. Therefore, in the Granada and Venezuela Basins, MHWs are less frequent but last longer. The seasonality of all the MHWs metrics in the spatial average are coherent with the first row of Figure 4.

The spatial assessment of the MHWs metrics was also performed using the CESM2 and 5MMM models to assess their ability to accurately reproduce the annual and seasonal spatial behavior (Figures S1 and S2, Supplementary Material. Other models' results shown in Figures S3–S6). When compared to observations (Figure 5), both models, namely CESM2 and 5MMM, underestimated the mean annual and seasonal intensity, with nodes' values ranging from 0.47 to 0.92 °C in the former and 0.34 to 0.58 °C in the latter. Furthermore, seasonal differences are smaller in CESM2 but nearly inexistant in 5MMM, which overall shows a poorer performance in this metric. However, a larger mean intensity can be observed in the Guajira upwelling region in both models, as also seen in observations.

Both models also underestimate the average number of events in the annual and seasonal frequency, ranging from 0.41 to 2.27 events per month in the CESM2 model and 0.42 to 2.67 events per month in the 5MMM model. However, both models correctly capture the larger number of events per month during the rainy season when compared to the dry season. Overall, 5MMM shows better performance in this metric.

The least accurate performance is in the MHWs duration, as both CESM2 and 5MMM models largely overestimate events' duration mean in all three cases. The MHWs mean duration ranges from 14.12 to 68.33 days in the CESM2 model and 10.01 to 83.33 days in the 5MMM model when the annual and two seasons' results were assessed. Moreover, seasonal and spatial differences do not follow observations (Figure 5, third row), as the models do not reproduce the weak spatial pattern in the MHWs frequency and duration, where MHWs in the eastern Caribbean are less frequent but with larger duration.

### 3.2.4. MHWs Events

We further investigated MHWs months to assess individual events in the CAR based on the NOAA OISST V2 dataset. We selected two extreme months from the baseline period (1983–2012). First, we assessed the month with the highest mean intensity (2.18 °C in February 1985) (Figure 4D). In this month, the largest mean SST anomaly occurred in the node located at 12.5 °N and 69.5 °W (2.32 °C). Only a neighbor node also recorded MHWs in this month but with a lower mean SST anomaly (2.03 °C). This result indicates that MHWs can evolve and impact very small areas, reaching high SST anomalies (Figure 6A).

The 1985 daily time series indicates that this extreme MHW had a 6-day duration (SST above the daily 90th percentile threshold—Figure 6B). Note that in this node, three short MHWs occurred in 1985. The other two MHW occurred at the end of April and in June, with a mean intensity (duration) of 2.38 °C (14 days) and 1.83 °C (13 days), respectively. Due to the methodology we used, MHWs do not include the seasonal SST variation, which can be clearly seen in the threshold time series. We did not include seasonality, as we want to better understand the events that force MHWs in the CAR. Nonetheless, due to the seasonal SST variations, less intense MHWs such as the one in June will have a larger impact, e.g., in the biosphere, than the more intense MHW that occurred at the beginning of February, as it finally reached the highest SST in 1985 for that node (30.17 °C). Therefore, MHWs that occur in the rainy season (SON), which seems to be more intense and frequent (Figure 5), will have more devastating impacts on the biosphere (e.g., [10]) but also would affect the hurricane season that peaks in this time of the year [63].
<figure><img src="image_10.png" /><figcaption>**Figure 6.** Extreme MHWs months in the Caribbean Sea. First row shows February 1985, the month with maximum mean intensity, while the second row shows March 2010, the month with largest mean duration. First column shows the mean (A) intensity and (C) duration from the nodes in which MHWs occur during the month. The white cross indicates the node whose daily SST time series (black line) is shown in (B,D) for the entire year to highlight seasonal variations. The time series includes the node's 30-year climatology (median) and 90th percentile, which is the threshold used to define MHWs. SST anomalies of the MHWs during the assessed months are highlighted in red, while other MHWs are visually represented in orange.</figcaption></figure>

Secondly, we selected the month with the longest mean duration (20.72 days in March 2010) (Figure 4F). In contrast, the results show MHWs occurring in nearly all nodes in the eastern Caribbean, lasting between 10 and 31 days (Figure 6C), which is coherent with the annual duration shown in Figure 5. MHWs with a larger duration occur in the north-eastern Caribbean, which coincides with an area of severe coral bleaching reported in the warm 2005 [10].

The 2010 daily time series from one of the nodes with the longest duration (13.5 °N, 69.5 °W) shows that the same MHW lasted 50 days, from February 20 to April 11 (Figure 6D). This allows us to explain a limitation in the methodology used in Figure 4, where we report the MHWs monthly mean duration, which does not account for MHWs that occur in more than one month. In this node and year, eight MHWs occurred with a duration (mean intensity) between 11 to 51 days (0.81 °C to 1.10 °C). Note again the importance of seasonality in the impacts, as the MHW in September, with a small intensity, reached the highest SST in 2010 for that node (30.18 °C). Additionally, it seems that seasonality also affects MHWs duration, as in the colder first quarter of the year, SST anomalies above the threshold last for more days. This behavior is also seen in Figures 4C and 5 in the third row. Overall, 2010 was a warm year, as most SST time series were above climatology in the assessed node.

A comparison of the MHWs coverage and other metrics in February 1985 and March 2010 shows how differently these events can evolve and impact the CAR. In February 1985, an intense MHW locally impacted a small region north of the Venezuela coast. Also, during this year, three short but intense MHWs affected this node. In contrast, many MHWs affected most of the eastern Caribbean with a mild intensity in March 2010. Note how, in this year, duration decreased from east to west, which would indicate a forcing entering the Caribbean from the northern end of the Lesser Antilles toward the west, producing the MHWs. This basin-wide warming, which was observed through the first half of the
year, seems to be linked to the tropical North Atlantic. We speculate that MHWs with large coverage in the Caribbean Sea are probably linked to seasonal and interannual variations of the Atlantic Warm Pool [65], whose intensity and area variations are linked to El Niño Southern Oscillation and North Atlantic Oscillation [66].

A reviewer suggested that the February 1985 MHW event might be spurious, as it was registered in only two nodes. This could be probably related to an OISST dataset limitation, as it is known that SST satellite products are less reliable close to the coasts [67]. To assess this suggestion, we analyzed the 10 most intense MHWs events. Results revealed that these intense events typically persist less than 12 days, with temperature anomalies >2 °C, with a spatial concentration within the 11.5–12.5°N and 69.5–73.5°W coordinates, in a region of active coastal upwelling. Some of these events appear to cover a small number of nodes, including an event in June 2011, when better SST data were available. Therefore, we conclude that the low-coverage February 1985 MHW event is probably a real event registered by the OISST dataset.

### 3.3. MHWs Projections for the 21st Century

CMIP6 model simulations were used to analyze expected climate change impacts to the CAR MHWs from 2040 to 2100, differentiating the range into two 30-year periods representing mid-century (2041–2070) and end-century (2071–2100) conditions. Trends in the MHWs metrics (mean intensity, frequency, and mean duration) were assessed under two SSPs scenarios, following the same methodology used to assess the 1983–2012 reference period. In this assessment, we prioritized the description of the 5MMM results because that dataset better reproduces the observed SST behavior when compared to CESM2, with lower normalized standard deviation and RMSE but higher correlation (Table 2, Figure 2). Furthermore, 5MMM better reproduces the OISST spatial and seasonal variability when compared to CESM2 (Section 3.2.1, Figure 3) (Section 3.2.3). However, results from the CESM2 are also included in Table 4.

<table><caption>**Table 4.** MHWs properties for mid-century (2041–2070) and end-century (2071–2100) conditions in the Caribbean Sea from the 5MMM and CESM2 models annual time series. Intensity and duration from spatially averaged annual values. Frequency from the sum of events per year from all Caribbean nodes. Metrics include the 30-year mean, standard deviation (Std), and the yearly maximum value (MaV).</caption><thead><tr><th rowspan="3"></th><th colspan="12">SSP245</th></tr><tr><th colspan="6">2041–2070</th><th colspan="6">2071–2100</th></tr><tr><th colspan="2">Intensity (°C)</th><th colspan="2">Frequency (Total Events)</th><th colspan="2">Duration (Days)</th><th colspan="2">Intensity (°C)</th><th colspan="2">Frequency (Total Events)</th><th colspan="2">Duration (Days)</th></tr></thead><tbody><tr><th>Mean</th><td>0.72</td><td>0.78</td><td>417.93</td><td>291.53</td><td>85.62</td><td>90.48</td><td>1.18</td><td>1.34</td><td>237.35</td><td>200.17</td><td>101.92</td><td>110.58</td></tr><tr><th>Std</th><td>0.18</td><td>0.34</td><td>105.82</td><td>122.40</td><td>10.61</td><td>19.80</td><td>0.12</td><td>0.32</td><td>25.54</td><td>12.20</td><td>7.54</td><td>4.60</td></tr><tr><th>MaV</th><td>1.05</td><td>1.71</td><td>644.20</td><td>625.00</td><td>103.64</td><td>109.81</td><td>1.40</td><td>1.89</td><td>302.40</td><td>238.00</td><td>117.30</td><td>115.00</td></tr></tbody><tbody><tr><th colspan="12">SSP585</th></tr><tr><th>Mean</th><td>1.10</td><td>1.29</td><td>299.47</td><td>239.40</td><td>95.66</td><td>98.30</td><td>2.23</td><td>2.75</td><td>196.71</td><td>186.79</td><td>111.84</td><td>130.56</td></tr><tr><th>Std</th><td>0.31</td><td>0.51</td><td>94.85</td><td>100.47</td><td>12.42</td><td>16.51</td><td>0.34</td><td>0.48</td><td>2.82</td><td>16.79</td><td>5.94</td><td>4.57</td></tr><tr><th>MaV</th><td>1.53</td><td>2.14</td><td>510.49</td><td>604.00</td><td>116.00</td><td>115.33</td><td>2.69</td><td>3.56</td><td>204.60</td><td>198.00</td><td>118.33</td><td>154.22</td></tr></tbody></table>

The MHWs mean intensity is projected to increase 0.146 (0.368) °C per decade in the 2041–2100 period under the SSP245 (SSP585) scenario (Figure 7A). Trends are larger for the end-century period (0.370 °C per decade) than for the mid-century period (0.322 °C per decade) under the SSP585 scenario. Conversely, mean intensity trends are larger for the mid-century period (0.152 °C per decade) than for the end-century period (0.102 °C per decade).
per decade) under SSP245, as this scenario reduces atmospheric greenhouse gases in the second half of the 21st century [68].

<figure><img src="image_11.png" /><figcaption>**Figure 7.** Caribbean Sea annual time series of MHWs metrics under SSP245 (blue) and SSP585 (red) socioeconomic pathways for the 5MMM ensemble (solid thick lines): (A) spatially averaged mean intensity, (B) frequency (sum of events per year from all Caribbean nodes), and (C) spatially averaged mean duration. Dashed lines represent the 95% confidence interval; 1983–2012 is the baseline period used to calculate the daily SST threshold in each node and used to define MHWs.</figcaption></figure>

The MHWs frequency shows negative trends regardless of the SSP scenario considered (Figure 7B). In the 2040–2100 period, a greater reduction is projected under SSP245 (–57.79 events per decade) when compared to SSP585 scenario (–37.39 events per decade). In both cases, at the end of the century, ~200 total events per year in the entire Caribbean are expected (~1 annual event per year for each of the 195 nodes). This reduction in the MHWs frequency is related to the projected increase in the events’ duration. In the 2040–2100 period, positive trends of 5.0 and 5.7 days per year are expected in the mean events’ duration under the SSP245 and SSP585 scenarios, respectively (Figure 7C). These projections indicate that for 2100, the mean MHWs duration in the CAR would be ~100 days. Clearly, the projected MHWs frequency reduces, as the events’ projected duration would stand for over three months. Note that in both metrics, the 95% confidence interval reduces toward the end of the century.

These results must be treated with caution for several reasons. First, the MHWs projected long duration is partially because the 90th percentile threshold used to define
MHWs is based on the SST conditions during the 1983–2012 period. Second, when we compared the model and observed the MHWs time series in the 1983–2012 period, we observed that the models underestimated the mean intensity and frequency, while they overestimated the mean duration (Figure 3). Therefore, projected MHWs metrics using climate models should not be interpreted in absolute terms; still, trends are useful to indicate expected future conditions in the CAR related to these extreme events. Following this line of thought, the results indicate that MHWs will be nearly permanent and more intense by the end of the century in the study area, based on the observed SST conditions from the 1983–2012 baseline period.

We also compared MHWs metrics computed for mid-century and end-century conditions, using both the CESM2 and 5MMM models' projections in the SSP245 and SSP585 scenarios (Table 4). These results can be compared to the MHWs metrics observed for the 1983–2012 period (Table 3).

In both scenarios, the mean intensity is higher in the 2071–2100 period when compared to the 2041–2070 period. More intense events are projected in the SSP585 scenario, with mean intensity values up to 2.2 and 2.8 °C in the end-century conditions for the 5MMM and CESM2 models, respectively. In all cases, CESM2 projects more intense MHWs than the 5MMM. Note that some years can reach mean intensity values of 2.7 and 3.6 °C in both models (SSP585). This indicates that the projected mean MHWs intensity in the CAR can be over twice the observed mean intensity in the 1983–2012 period (0.94 °C in Table 3), with years that might reach SST anomalies >3 °C with respect to the reference period.

In all cases, mean duration is larger in the SSP585 scenario when compared to the SSP245 as well as in the end-century when compared to mid-century period (Table 4). Conversely, annual dispersion (standard deviation) is smaller in all cases in the end-century than in the mid-century conditions. CESM2 in all cases projects a larger mean MHWs duration when compared to the 5MMM ensemble. In the 1983–2012 reference period, both models nearly doubled the observed MHWs duration (Table 3); still, the models project that for 2071–2100, the mean MHWs duration in CAR will be over six times the modeled mean duration in the reference period.

As a consequence of the large increase in the mean MHWs duration, the frequency will decrease with time, as also seen in the extreme radiative scenario when compared to the SSP245 pathway (Table 4). Similarly, as CESM2 projects larger mean durations, it also projects less frequent events when compared to 5MMM.

Overall, the comparison between the 5MMM and CESM2 models shows differences in their MHWs metrics projections, with CESM2 indicating more intense, less frequent, and longer-lasting events, particularly under the SSP585 scenario.

## 4. Summary and Final Remarks

In this study, we assessed MHWs in the CAR during the 1983–2012 period, using the NOAA OISST V2 daily SST dataset. MHWs were assessed based on three metrics, namely mean intensity, mean duration, and its frequency, which were used to describe their interannual, annual, and sub annual variability as well as their spatial behavior from the annual mean and the dry and rainy seasons. In addition, cases from two extreme months (higher mean intensity and longest mean duration) were included to show how different these events can be in the region. In this study, SST seasonality did not contribute to MHW events.

To evaluate future MHW events in the Caribbean Sea and their relationship with global warming, we used results from 19 CMIP6 models during the historical period (1983–2012) and 2014–2100 projections under two socioeconomic pathways (SSP245 and SSP585). Daily SST time series from the models’ historical period were compared to observations using a Taylor diagram (Figure 2). The best results were obtained from a multi-model mean using the five best-performing models and by the CESM2 model. Therefore, these results are used through the paper.
In the 1983–2012 period, mean duration and frequency from the observed MHWs showed significant trends but not in their intensity. The latter was unexpected, as SST trends in the CAR were reported for the 1960–2005 and 2005–2050 periods [60]. However, this behavior seems to be related to a small difference in the probability density function of extreme SST values (those included in the MHWs events) when comparing the first and second 15-year period of the 1983–2012 time span. In the same period as the observations, the models underestimated intensity and frequency while overestimating duration; however, they can show also significant trends in mean duration and frequency (Figure 3). In addition, the models showed limitations in reproducing the seasonality of the MHWs metrics (Figures S1 and S2, Supplementary Material).

Some studies have shown that the coarse resolution of GCMs leads to longer-lasting MHWs due an underestimation of SST anomalies related to the ability to reproduce small-scale features, producing spatially and temporally smoothed SST fields [6,21]. This weakness is preserved in the CMIP6 models despite the improvements in resolution.

Observations in the 1983–2012 period also showed large interannual and monthly variability in the MHWs metrics (Figure 4). For instance, 2010 was a year with many long-lasting MHWs, while years such as 1991 had nine consecutive months with no events. Interannual variability was largest in the MHWs frequency (standard deviation in Table 3).

Observed MHWs seasonality was weak and different in the three used metrics. The rainy season showed more intense and frequent MHWs, while the dry (cold and windy) season had longer-lasting MHWs, which corresponds to the inverse relationship between frequency and duration (Figures 4 and 5). As was also supported by two cases (Figure 6), MHWs in the CAR can be very different in their spatial coverage as well as in their mean intensity and duration. With these cases, it was shown that more intense and frequent MHWs in the rainy season will produce larger impacts, as this period coincides with the warmer months in the year.

Limited studies have been conducted about MHWs in the Caribbean Sea [11,32] using the same methodology we used, reporting mean duration trends of 0.20 days per decade, while our results showed a higher positive significant trend of $1.47 \pm 0.29$ days per decade. However, ref. [32] considered a longer time frame (1981–2018), evaluating specific coral reef locations in the CAR. Although the other MHWs metrics we used are not directly comparable to the results from [32], they also concluded that MHWs in recent decades have increased in their occurrence and duration (see Figure 3) [6,32]. Furthermore, the spatially averaged MHWs values from 1983 to 2012 obtained in our study ($1.00 \pm 0.12$ °C for mean intensity, $1.86 \pm 0.23$ events per month in frequency, and $11.31 \pm 2.11$ days for mean duration) are in the limits of the reported values in the CAR from a global assessment performed by [11].

In the CAR, significant positive trends are projected for mean duration and intensity and negative trends for frequency, regardless of the radiative scenario used; however, the trends are larger in the SSP585 scenario (Figure 7). The results obtained are similar to other global ocean studies relying on observations and CMIP5/CMIP6 GCMs [8,12,21]. Note that due to the increase in intensity with time, by the end of the 21st century, MHWs could reach anomalies >3 °C, with respect to the threshold (90th percentile) established in the 1983–2012 period (Table 4).

Models' projections indicate that by the end of the 21st century, actual extreme MHWs conditions in the CAR could be a regular condition, as the mean event duration will largely increase (frequency reduction). Moreover, mean intensity will be over twice the intensity observed in the 1983–2012 reference period. These conditions will probably have strong implications for the ocean biosphere and affect the hurricane season, as already seen during warm years (e.g., [10,63]). Nonetheless, MHWs projections are sensitive to the radiative scenario, and as seen in the SSP245, MHWs are expected to be less intense and with shorter duration than in the SSP585 pathway. Therefore, efforts to reduce atmospheric greenhouse gases should continue and strengthen.

Thus, the main important conclusions can be summarized in these points:
- CMIP6 models showed significant trends in the spatially averaged MHWs frequency and mean duration but not in mean intensity during 1983–2012;

- MHWs are more intense and frequent in warm–rainy months and longer-lasting in late boreal winter;

- MHWs conditions at the beginning of the century will be nearly permanent in the Caribbean’s future.

Finally, these results highlight the importance of continued study of MHWs in the CAR since in this region, there are a vast number of coral reefs, which are home to about 25% of all marine fish species [69], and as such, they are an important spotlight for fishing and coastal protection [70]. These species are very sensitive to changes in temperature, light, ultraviolet radiation, salinity, pH, and excessive sedimentation, being strongly affected by changes in the thermal environment and other stressors, including fishing and pollution [70]. The projected increase in duration and intensity of MHWs could cause a dramatic shift in the composition and functioning of Caribbean coral reef ecosystems, putting them in danger of coral bleaching, with a high risk of disappearance [32]. However, there are several other factors to consider regarding the impact of climate change on coral reefs, such as inter-annual oscillations like the ENSO and interactions between ecological communities, among others [71–73]. How have coral reefs been impacted, and how can MHWs affect these ecosystems in the CAR? In what manner do MHWs affect coral reefs in the CAR under different protection regimes? How are other factors, such as anthropogenic activity, affecting or not affecting the reefs? These are questions that are not yet completely answered.

Understanding how MHWs can evolve in the future in the CAR holds significance for effective marine resource management, necessitating a thorough comprehension of their drivers. Strategic decisions in resource management often rely on climate forecasts, with ecosystem-based approaches requiring consideration of environmental stressors [74]. The results from this study can be used as a starting point for understanding these questions and help stakeholders to make good decisions to protect the Caribbean environment in the future.

**Supplementary Materials:** The following supporting information can be downloaded at: https://www.mdpi.com/article/10.3390/rs16132357/s1, Figure S1: MHWs spatial behavior from the CESM2 model in the CAR; Figure S2: MHWs spatial behavior from the 5MMM model in the CAR; Figure S3: MHWs spatial behavior from the ACCESS-ESM1-5 model in the CAR; Figure S4: MHWs spatial behavior from the GFDL-ESM4 model in the CAR; Figure S5: MHWs spatial behavior from the NorESM2-LM model in the CAR; Figure S6: MHWs spatial behavior from the NorESM2-MM model in the CAR; Table S1: Seasonal cycle assessment; Table S2: MHWs properties for the 1983–2012 period in the Caribbean Sea from OISSTV2, CESM2, ACCESS-ESM1-5, GFDL-ESM4, NorESM2-LM and NorESM2-MM annual time series.

**Author Contributions:** Conceptualization, D.F.B.U. and R.R.T.P.; methodology, D.F.B.U. and L.R.-L.; software, D.F.B.U.; validation, D.F.B.U., R.R.T.P. and L.R.-L.; formal analysis, D.F.B.U.; investigation, D.F.B.U.; resources, L.B.; data curation, D.F.B.U.; writing—original draft preparation, D.F.B.U.; writing—review and editing, R.R.T.P., L.B, L.R.-L. and M.C.A.; visualization, D.F.B.U.; supervision, R.R.T.P. and L.R.-L.; project administration, R.R.T.P.; funding acquisition, L.B. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research received no external funding.

**Data Availability Statement:** Available upon request from the corresponding authors.

**Acknowledgments:** We acknowledge the Coupled Model Intercomparison Project (CMIP) launched by the World Climate Research Program (WCRP) for CMIP6 data.

**Conflicts of Interest:** The authors declare no conflicts of interest.
## References

1. Chen, W.; Jiang, Z.; Li, L. Probabilistic Projections of Climate Change over China under the SRES A1B Scenario Using 28 AOGCMs. *J. Clim.* **2011**, *24*, 4741–4756. [CrossRef]

2. Roemmich, D.; John Gould, W.; Gilson, J. 135 Years of Global Ocean Warming between the Challenger Expedition and the Argo Programme. *Nat. Clim. Chang.* **2012**, *2*, 425–428. [CrossRef]

3. Cheng, L.; Zhu, J.; Abraham, J.; Trenberth, K.E.; Fasullo, J.T.; Zhang, B.; Yu, F.; Wan, L.; Chen, X.; Song, X. 2018 *Continues Record Global Ocean Warming*; Springer: Berlin/Heidelberg, Germany, 2019; ISBN 1861-9533.

4. Johnson, G.C.; Lyman, J.M. Warming Trends Increasingly Dominate Global Ocean. *Nat. Clim. Chang.* **2020**, *10*, 757–761. [CrossRef]

5. Hobday, A.J.; Alexander, L.V.; Perkins, S.E.; Smale, D.A.; Straub, S.C.; Oliver, E.C.; Benthuysen, J.A.; Burrows, M.T.; Donat, M.G.; Feng, M. A Hierarchical Approach to Defining Marine Heatwaves. *Prog. Oceanogr.* **2016**, *141*, 227–238. [CrossRef]

6. Oliver, E.C.; Donat, M.G.; Burrows, M.T.; Moore, P.J.; Smale, D.A.; Alexander, L.V.; Benthuysen, J.A.; Feng, M.; Sen Gupta, A.; Hobday, A.J. Longer and More Frequent Marine Heatwaves over the Past Century. *Nat. Commun.* **2018**, *9*, 1–12. [CrossRef] [PubMed]

7. Frölicher, T.L.; Laufkötter, C. Emerging Risks from Marine Heat Waves. *Nat. Commun.* **2018**, *9*, 1–4. [CrossRef] [PubMed]

8. Smale, D.A.; Wernberg, T.; Oliver, E.C.; Thomsen, M.; Harvey, B.P.; Straub, S.C.; Burrows, M.T.; Alexander, L.V.; Benthuysen, J.A.; Donat, M.G. Marine Heatwaves Threaten Global Biodiversity and the Provision of Ecosystem Services. *Nat. Clim. Chang.* **2019**, *9*, 306–312. [CrossRef]

9. Pearce, A.; Jackson, G.; Moore, J.; Feng, M.; Gaughan, D.J. *The "Marine Heat Wave" off Western Australia during the Summer of 2010/11*; Western Australian Fisheries and Marine Research Laboratories: Hillarys, WA, Australia, 2011.

10. Eakin, C.M.; Morgan, J.A.; Heron, S.F.; Smith, T.B.; Liu, G.; Alvarez-Filip, L.; Baca, B.; Bartels, E.; Bastidas, C.; Bouchon, C. Caribbean Corals in Crisis: Record Thermal Stress, Bleaching, and Mortality in 2005. *PLoS ONE* **2010**, *5*, e13969. [CrossRef] [PubMed]

11. Holbrook, N.J.; Scannell, H.A.; Sen Gupta, A.; Benthuysen, J.A.; Feng, M.; Oliver, E.C.; Alexander, L.V.; Burrows, M.T.; Donat, M.G.; Hobday, A.J. A Global Assessment of Marine Heatwaves and Their Drivers. *Nat. Commun.* **2019**, *10*, 1–13. [CrossRef]

12. Plecha, S.M.; Soares, P.M. Global Marine Heatwave Events Using the New CMIP6 Multi-Model Ensemble: From Shortcomings in Present Climate to Future Projections. *Environ. Res. Lett.* **2020**, *15*, 124058. [CrossRef]

13. IPCC The Physical Science Basis. *Contribution of Working Group I to the Sixth Assessment Report of the Intergovernmental Panel on Climate Change 2021*; IPCC: Geneva, Switzerland, 2021.

14. McCarthy, J.J. *Climate Change 2001: Impacts, Adaptation, and Vulnerability: Contribution of Working Group II to the Third Assessment Report of the Intergovernmental Panel on Climate Change*; Cambridge University Press: Cambridge, UK, 2001; Volume 2, ISBN 0-521-01500-6.

15. O'Neill, B.C.; Tebaldi, C.; van Vuuren, D.P.; Eyring, V.; Friedlingstein, P.; Hurtt, G.; Knutti, R.; Kriegler, E.; Lamarque, J.-F.; Lowe, J.; et al. The Scenario Model Intercomparison Project (ScenarioMIP) for CMIP6. *Geosci. Model Dev.* **2016**, *9*, 3461–3482. [CrossRef]

16. Eyring, V.; Bony, S.; Meehl, G.; Senior, C.; Stevens, B.; Ronald, S. Overview of the Coupled Model Intercomparison Project Phase 6 (CMIP6) Experimental Design and Organisation. *Geosci. Model Dev. Discuss.* **2015**, *8*, 10539–10583. [CrossRef]

17. Qiu, Z.; Qiao, F.; Jang, C.J.; Zhang, L.; Song, Z. Evaluation and Projection of Global Marine Heatwaves Based on CMIP6 Models. *Deep Sea Res. Part II Top. Stud. Oceanogr.* **2021**, *194*, 104998. [CrossRef]

18. Pilo, G.S.; Holbrook, N.J.; Kiss, A.E.; Hogg, A.M. Sensitivity of Marine Heatwave Metrics to Ocean Model Resolution. *Geophys. Res. Lett.* **2019**, *46*, 14604–14612. [CrossRef]

19. Richardson, P.L. Caribbean Current and Eddies as Observed by Surface Drifters. *Deep Sea Res. Part II Top. Stud. Oceanogr.* **2005**, *52*, 429–463. [CrossRef]

20. Darmaraki, S.; Somot, S.; Sevault, F.; Nabat, P.; Cabos Narvaez, W.D.; Cavicchia, L.; Djurdjevic, V.; Li, L.; Sannino, G.; Sein, D.V. Future Evolution of Marine Heatwaves in the Mediterranean Sea. *Clim. Dyn.* **2019**, *53*, 1371–1392. [CrossRef]

21. Oliver, E.C.; Burrows, M.T.; Donat, M.G.; Sen Gupta, A.; Alexander, L.V.; Perkins-Kirkpatrick, S.E.; Benthuysen, J.A.; Hobday, A.J.; Holbrook, N.J.; Moore, P.J. Projected Marine Heatwaves in the 21st Century and the Potential for Ecological Impact. *Front. Mar. Sci.* **2019**, *6*, 734. [CrossRef]

22. Plecha, S.M.; Soares, P.M.; Silva-Fernandes, S.M.; Cabos, W. On the Uncertainty of Future Projections of Marine Heatwave Events in the North Atlantic Ocean. *Clim. Dyn.* **2021**, *56*, 2027–2056. [CrossRef]

23. Taylor, K.E.; Stouffer, R.J.; Meehl, G.A. An Overview of CMIP5 and the Experiment Design. *Bull. Am. Meteorol. Soc.* **2012**, *93*, 485–498. [CrossRef]

24. Chen, Z.; Zhou, T.; Zhang, L.; Chen, X.; Zhang, W.; Jiang, J. Global Land Monsoon Precipitation Changes in CMIP6 Projections. *Geophys. Res. Lett.* **2020**, *47*, e2019GL086902. [CrossRef]

25. Grose, M.R.; Narsey, S.; Delage, F.P.; Dowdy, A.J.; Bador, M.; Boschat, G.; Chung, C.; Kajtar, J.B.; Rauniyar, S.; Freund, M.B. Insights from CMIP6 for Australia's Future Climate. *Earth's Future* **2020**, *8*, e2019EF001469. [CrossRef]

26. Kajtar, J.B.; Hernaman, V.; Holbrook, N.J.; Petrelli, P. Tropical Western and Central Pacific Marine Heatwave Data Calculated from Gridded Sea Surface Temperature Observations and CMIP6. *Data Brief* **2022**, *40*, 107694. [CrossRef] [PubMed]

27. Liu, H.; Song, Z.; Wang, X.; Misra, V. *An Ocean Perspective on CMIP6 Climate Model Evaluations*; Elsevier: Amsterdam, The Netherlands, 2022; p. 105120. ISBN 0967-0645.
28. Riahi, K.; van Vuuren, D.P.; Kriegler, E.; Edmonds, J.; O'Neill, B.C.; Fujimori, S.; Bauer, N.; Calvin, K.; Dellink, R.; Fricko, O.; et al. The Shared Socioeconomic Pathways and Their Energy, Land Use, and Greenhouse Gas Emissions Implications: An Overview. *Glob. Environ. Chang.* **2017**, *42*, 153–168. [CrossRef]

29. Samuels, T.; Rynearson, T.A.; Collins, S. Surviving Heatwaves: Thermal Experience Predicts Life and Death in a Southern Ocean Diatom. *Front. Mar. Sci.* **2021**, *8*, 600343. [CrossRef]

30. Holbrook, N.J.; Hernaman, V.; Koshiba, S.; Lako, J.; Kajtar, J.B.; Amosa, P.; Singh, A. Impacts of Marine Heatwaves on Tropical Western and Central Pacific Island Nations and Their Communities. *Glob. Planet. Chang.* **2022**, *208*, 103680. [CrossRef]

31. Donovan, M.K.; Burkepile, D.E.; Kratochwill, C.; Shlesinger, T.; Sully, S.; Oliver, T.A.; Hodgson, G.; Freiwald, J.; van Woesik, R. Local Conditions Magnify Coral Loss after Marine Heatwaves. *Science* **2021**, *372*, 977–980. [CrossRef]

32. Bove, C.B.; Mudge, L.; Bruno, J.F. A Century of Warming on Caribbean Reefs. *PLoS Clim. Change*, **2022**, *1*, e000002. [CrossRef]

33. Smith, K.E.; Burrows, M.T.; Hobday, A.J.; Sen Gupta, A.; Moore, P.J.; Thomsen, M.; Wernberg, T.; Smale, D.A. Socioeconomic Impacts of Marine Heatwaves: Global Issues and Opportunities. *Science* **2021**, *374*, eabj3593. [CrossRef]

34. Reyes, O.; Manta, G.; Carrillo, L. Marine Heatwaves and Marine Cold-Spells on the Yucatan Shelf-Break Upwelling Region. *Cont. Shelf Res.* **2022**, *239*, 104707.

35. Martinez, S. THEMATIC REPORT FOR THE CENTRAL/SOUTH AMERICAN SUB-REGION. In *Unpublished. CLME Project Implementation Unit, Centre for Resource Management and Environmental Studies (CERMES)*; University of the West Indies: Saint Michael, Barbados, 2007.

36. Wilkinson, T.A.; Wiken, E.; Creel, J.B.; Hourigan, T.F.; Agardy, T. Marine Ecoregions of North America; Instituto Nacional de Ecología: Ciudad de México, Mexico, 2009; ISBN 2-923358-41-4.

37. Torres, R.R.; Latandret, S.; Salon, J.; Dagua, C. Water Masses in the Caribbean Sea and Sub-Annual Variability in the Guajira Upwelling Region. *Ocean Dyn.* **2023**, *73*, 39–57. [CrossRef]

38. GEBCO Gridded Bathymetry Data. Gridded Bathymetry Data. Available online: https://www.gebco.net/data_and_products/gridded_bathymetry_data/ (accessed on 30 April 2021).

39. Jacox, M.G. Marine Heatwaves in a Changing Climate; Nature Publishing Group: London, UK, 2019.

40. Hayashida, H.; Matear, R.J.; Strutton, P.G.; Zhang, X. Insights into Projected Changes in Marine Heatwaves from a High-Resolution Ocean Circulation Model. *Nat. Commun.* **2020**, *11*, 4352. [CrossRef] [PubMed]

41. Sen, P. Estimates of the Regression Coefficient Based on Kendall's Tau. *J. Am. Stat. Assoc.* **1968**, *63*, 1379–1389. [CrossRef]

42. Kendall, M.G. Rank Correlation Methods. Griffin, London. *J. Econ.* **1975**, *13*, 245–259.

43. Mann, H.B. Nonparametric Tests against Trend. *Econom. J. Econom. Soc.* **1945**, *13*, 245–259. [CrossRef]

44. Mitchell, J.M.; Dzerdzeevskii, B.; Flohn, H. *Climate Change*, 79th ed.; World Metereological Organization: Geneva, Switzerland, 1966.

45. Huang, B.; Liu, C.; Banzon, V.; Freeman, E.; Graham, G.; Hankins, B.; Smith, T.; Zhang, H.-M. Improvements of the Daily Optimum Interpolation Sea Surface Temperature (DOISST) Version 2.1. *J. Clim.* **2021**, *34*, 2923–2939. [CrossRef]

46. Manta, G.; de Mello, S.; Trinchin, R.; Badagian, J.; Barreiro, M. The 2017 Record Marine Heatwave in the Southwestern Atlantic Shelf. *Geophys. Res. Lett.* **2018**, *45*, 12449–12456. [CrossRef]

47. Perkins-Kirkpatrick, S.; King, A.D.; Cougnon, E.A.; Grose, M.R.; Oliver, E.C.J.; Holbrook, N.; Lewis, S.; Pourasghar, F. The Role of Natural Variability and Anthropogenic Climate Change in the 2017/18 Tasman Sea Marine Heatwave. *Bull. Am. Meteorol. Soc.* **2019**, *100*, S105–S110. [CrossRef]

48. Ibrahim, O.; Mohamed, B.; Nagy, H. Spatial Variability and Trends of Marine Heat Waves in the Eastern Mediterranean Sea over 39 Years. *J. Mar. Sci. Eng.* **2021**, *9*, 643. [CrossRef]

49. Cui, T.; Zhang, J.; Groom, S.; Sun, L.; Smyth, T.; Sathyendranath, S. Validation of MERIS Ocean-Color Products in the Bohai Sea: A Case Study for Turbid Coastal Waters. *Remote Sens. Environ.* **2010**, *114*, 2326–2336. [CrossRef]

50. Taylor, K.E.; Juckes, M.; Balaji, V.; Cinquini, L.; Denvil, S.; Durack, P.J.; Elkington, M.; Guilyardi, E.; Kharin, S.; Lautenschlager, M. CMIP6 Global Attributes, DRS, Filenames, Directory Structure, and CV's; PCMDI Document; 2018. Available online: https://docs.google.com/document/d/1h0r8RZrf3-8egBMMh7aqLwy3snpD6_MrDz1q8n5XUk/edit (accessed on 26 June 2024).

51. Gidden, M.J.; Riahi, K.; Smith, S.J.; Fujimori, S.; Luderer, G.; Kriegler, E.; van Vuuren, D.P.; van den Berg, M.; Feng, L.; Klein, D.; et al. Global Emissions Pathways under Different Socioeconomic Scenarios for Use in CMIP6: A Dataset of Harmonized Emissions Trajectories through the End of the Century. *Geosci. Model Dev.* **2019**, *12*, 1443–1475. [CrossRef]

52. Torres, R.R.; Tsimplis, M.N. Seasonal Sea Level Cycle in the Caribbean Sea. *J. Geophys. Res. Ocean.* **2012**, *117*, C07011. [CrossRef]

53. Kim, Y. KIOST KIOST-ESM Model Output Prepared for CMIP6 CMIP Historical. *Earth Syst. Grid Fed. 2019*.

54. Taylor, K.E. Summarizing Multiple Aspects of Model Performance in a Single Diagram. *J. Geophys. Res. Atmos.* **2001**, *106*, 7183–7192. [CrossRef]

55. Wehner, M.; Gleckler, P.; Lee, J. Characterization of Long Period Return Values of Extreme Daily Temperature and Precipitation in the CMIP6 Models: Part 1, Model Evaluation. *Weather Clim. Extrem.* **2020**, *30*, 100283. [CrossRef]

56. Yang, X.; Huang, P. Improvements in the Relationship between Tropical Precipitation and Sea Surface Temperature from CMIP5 to CMIP6. *Clim. Dyn.* **2022**, *60*, 3319–3337. [CrossRef]

57. Nooni, I.K.; Ogou, F.K.; Chaibou, A.A.S.; Nakoty, F.M.; Gnitou, G.T.; Lu, J. Evaluating CMIP6 Historical Mean Precipitation over Africa and the Arabian Peninsula against Satellite-Based Observation. *Atmosphere* **2023**, *14*, 607. [CrossRef]

58. Gleckler, P.J.; Taylor, K.E.; Doutriaux, C. Performance Metrics for Climate Models. *J. Geophys. Res. Atmos.* **2008**, *113*. [CrossRef]
59. Flato, G. Evaluation of Climate Models. In *Climate Change 2013: The Physical Science Basis. Contribution of Working Group I to the Fifth Assessment Report of the Intergovernmental Panel on Climate Change*; IPCC: Geneva, Switzerland, 2013; pp. 741–866.

60. Bustos, D.; Torres, R. Ocean and Atmosphere Changes in the Caribbean Sea during the Twenty-First Century Using CMIP5 Models. *Ocean Dyn.* **2021**, *71*, 757–777. [CrossRef]

61. Costa, N.V.; Rodrigues, R.R. Future Summer Marine Heatwaves in the Western South Atlantic. *Geophys. Res. Lett.* **2021**, *48*, e2021GL094509. [CrossRef]

62. Le Grix, N.; Zscheischler, J.; Rodgers, K.B.; Yamaguchi, R.; Frölicher, T.L. Hotspots and Drivers of Compound Marine Heatwaves and Low Net Primary Production Extremes. *Biogeosciences* **2022**, *19*, 5807–5835. [CrossRef]

63. Trenberth, K.E.; Shea, D.J. Atlantic Hurricanes and Natural Variability in 2005. *Geophys. Res. Lett.* **2006**, *33*. [CrossRef]

64. Gamble, D.W.; Curtis, S. Caribbean Precipitation: Review, Model and Prospect. *Prog. Phys. Geogr. Earth Environ.* **2008**, *32*, 265–276. [CrossRef]

65. Wang, C.; Lee, S. Atlantic Warm Pool, Caribbean Low-Level Jet, and Their Potential Impact on Atlantic Hurricanes. *Geophys. Res. Lett.* **2007**, *34*. [CrossRef]

66. Enfield, D.B.; Lee, S.-K.; Wang, C. How Are Large Western Hemisphere Warm Pools Formed? *Prog. Oceanogr.* **2006**, *70*, 346–365. [CrossRef]

67. Boehme, L.; Lonergan, M.; Todd, C.D. Comparison of Gridded Sea Surface Temperature Datasets for Marine Ecosystem Studies. *Mar. Ecol. Prog. Ser.* **2014**, *516*, 7–22. [CrossRef]

68. O'Neill, B.C.; Kriegler, E.; Ebi, K.L.; Kemp-Benedict, E.; Riahi, K.; Rothman, D.S.; van Ruijven, B.J.; van Vuuren, D.P.; Birkmann, J.; Kok, K.; et al. The Roads Ahead: Narratives for Shared Socioeconomic Pathways Describing World Futures in the 21st Century. *Glob. Environ. Chang.* **2017**, *42*, 169–180. [CrossRef]

69. DellaSala, D.A.; Goldstein, M.I.; Elias, S.A.; Lacher, T.E.; Pyare, S.; Jennings, B.; Mineau, P. *Encyclopedia of the Anthropocene*; Elsevier: Amsterdam, The Netherlands, 2018; ISBN 0-12-813576-X.

70. Cramer, K.L.; Jackson, J.B.; Donovan, M.K.; Greenstein, B.J.; Korpanty, C.A.; Cook, G.M.; Pandolfi, J.M. Widespread Loss of Caribbean Acroporid Corals Was Underway before Coral Bleaching and Disease Outbreaks. *Sci. Adv.* **2020**, *6*, eaax9395. [CrossRef] [PubMed]

71. Graham, N.A.; Jennings, S.; MacNeil, M.A.; Mouillot, D.; Wilson, S.K. Predicting Climate-Driven Regime Shifts versus Rebound Potential in Coral Reefs. *Nature* **2015**, *518*, 94–97. [CrossRef]

72. Lough, J.M.; Anderson, K.D.; Hughes, T.P. Increasing Thermal Stress for Tropical Coral Reefs: 1871–2017. *Sci. Rep.* **2018**, *8*, 6079. [CrossRef]

73. Mumby, P.J.; Flower, J.; Chollett, I.; Box, S.J.; Bozec, Y.-M.; Fitzsimmons, C.; Forster, J.; Gill, D.; Griffith-Mumby, R.; Oxenford, H.A. *Towards Reef Resilience and Sustainable Livelihoods: A Handbook for Caribbean Coral Reef Managers*; University of Exeter: Exeter, UK, 2014; ISBN 0-902746-31-6.

74. Barbeaux, S.J.; Holsman, K.; Zador, S. Marine Heatwave Stress Test of Ecosystem-Based Fisheries Management in the Gulf of Alaska Pacific Cod Fishery. *Front. Mar. Sci.* **2020**, *7*, 703. [CrossRef]

**Disclaimer/Publisher's Note:** The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.