Marine Heatwaves in the South China Sea: Tempo-Spatial Pattern and Its Association with Large-Scale Circulation

Abstract

A marine heatwave (MHW) can significantly harm marine ecosystems and fisheries. Based on a remotely sensed sea surface temperature (SST) product, this study investigated MHWs behaviors in the South China Sea (SCS) throughout the warm season (May to September) from 1982 to 2020. The distributions of the three MHW indices used in this study showed significant latitudinal variations: more frequent, longer, and more intense MHWs appear in the northern SCS, and less frequent, shorter, and weaker MHWs appear in the southern SCS. Using the empirical orthogonal function (EOF) method, we found that the first leading modes of the three MHW indices account for more than half of the total variance. The first leading modes reveal uniform anomalies throughout the SCS, with the maximum in the deep central portion and its surroundings. Their corresponding time series showed significant interdecadal variations, with a turning point around 2009. Since 2010, the SCS has seen an increase in the frequency, length, and severity of MHWs. The incidence of MHWs has been linked to the presence of stable near-surface anticyclonic anomalies, which reduced cloud cover and increased solar radiation. This abnormal pattern was usually accompanied by the intensification and westward shift of the western North Pacific subtropical high (WNPSH). The findings imply that MHWs in the SCS may be predictable on interannual and decadal scales.

1. Introduction

The majority of the extra energy created by rising greenhouse gas (GHG) concentrations has been absorbed by the global oceans over the last century [1,2]. Studies showed that the additional energy warms the ocean, particularly the surface, resulting in marine high-temperature events known as marine heatwaves (MHWs) [3,4]. MHWs, in the last two decades, are becoming more frequent, longer, and more intense in most parts of the world [5,6,7]. More than 80% of the global oceans experienced at least one MHW in 2019 and 2020 [8]. MHWs have had a negative impact on marine organisms and ecosystems by pushing them to or beyond their thermal tolerance limits [4]. Research confirmed that long-distance migration of marine species, widespread coral bleaching, mass kelp forest mortality, and even economic losses in fishery industries have all been linked to MHWs [9,10]. MHWs are one of the primary indicators used in assessing the impact of climate change.

Recent studies have identified that the MHWs can be caused by the interaction of many local processes, large-scale climate modes, and teleconnections [4,11]. For example, atmospheric blocking over the ocean (e.g., Tasman Sea and Northwest Atlantic) can reduce cloud cover, enhance solar radiation, and calm winds, leading to MHWs [11]. In the context of global warming, atmospheric blockings have become more frequent and will be the key system in the future for the formation of large-scale MHWs [12]. Ocean dynamic processes, such as the enhancement of warm advection, shallowing of mixing layer, enhancement of ocean stratification, weakening of vertical mixing and upwelling, etc., are advantageous to the occurrence of MHWs [4]. In addition, MHWs can also be caused by El Niño and Southern Oscillation (ENSO), the Indian Ocean basin mode (IOB), and the Pacific Decadal Oscillation (PDO) in Indo-Pacific regions [13,14] remotely through atmospheric or oceanic teleconnections. In short, the MHWs in various regions can be controlled or modulated by a variety of drivers at various spatial and temporal scales [11].

As part of the Indo-Pacific warm pool, the South China Sea (SCS) (Figure 1a,b) is one of the world’s largest marginal seas, and it plays an important role in the economies and development of the surrounding nations due to its high marine biological productivity and important ecological functions [15]. Furthermore, The SCS rivals the Coral Triangle for the importance and diversity of its reefs. The SCS is home to an estimated 571 known species of reef corals [16]. The total reef area was estimated at about 37,935 km², accounting for about 5% of the world’s total coral reef area [17]. Corals in the SCS experienced unprecedented mass bleaching in 2015, resulting in 40% coral mortality [18]. In August-September of 2020, more widespread coral bleaching events occurred in the northern SCS, coastal areas of Guangdong, and Guangxi Provinces of China (Figure 1b), with up to 100% bleaching and an 86% mortality rate [19]. Coral bleaching has been linked to unusually severe MHWs in these areas. The research discovered that the marginal seas warmed much faster than the open oceans in a warming climate [20], indicating an increased likelihood of MHWs in these regions. Therefore, it is critical to investigate the characteristics, long-term variability, and drivers of MHWs in this region to improve prediction and build an early warning system.

Compared with other marginal seas, the SCS has fewer studies on MHWs [21]. The SCS has a typically tropical climate with long-lasting hot weather, and during the warm season, MHWs may easily exceed the thermal tolerance of marine organisms, leading to severe mass mortality [22]. Thus, we only considered MHWs that occurred during the 5-month warm season (May to September) in this region (more details are shown in Section 2). The spatiotemporal characteristics of warm season MHWs in the SCS are investigated for the period 1982–2020. The frequency, duration, and intensity of MHWs are the three key indicators analyzed. Furthermore, Empirical Orthogonal Function (EOF) is used to examine the dominant patterns of MHW and their regional characteristics. Here, we focus on the MHWs’ association with the large-scale atmospheric conditions which may initiate and sustain regional MHWs.

This paper is organized as follows. Section 2 describes the study region, data, and methods used. Section 3 depicts the results of our work, including climatology and long-term trends in annual MHWs indices, dominant Modes of MHWs variations, and connection with large-scale atmospheric circulation. Section 4 and Section 5 provide the discussion and conclusion.

2. Data and Methods

2.1. Study Region

The SCS is located in the tropical-to-subtropical western North Pacific, with geographical coordinates of 5°S~25°N and 100°~125°E (Figure 1a,b). It is bounded on the north by the shores of South China; on the east by the islands of Taiwan and northwestern Philippines; on the south by Borneo, eastern Sumatra, and the Bangka Belitung Islands. It covers about 3.6 × 106 km² and has a mean depth of 1212 m. The deep basin, known as the China Sea Basin, has a maximum depth of 5016 m in its northern central part. The SCS serves as a link between the East Asian and South Asian monsoon systems. The monsoon-dominated seasonal patterns distinguish it from other low-latitude waters that are insensitive to seasonal cycles; such characteristics drew the interest of researchers in its oceanic dynamic processes and long-term climate changes [23].

2.2. Data Description

Global high-resolution daily SST estimated by optimum interpolation (OI) of in situ and Advanced Very High-Resolution Radiometer infrared satellite SST data at 1/4° (OISST V2) was used in this study (the data can be downloaded from https://psl.noaa.gov/data/gridded/data.noaa.oisst.v2.highres.html, accessed on 10 April 2021) [24]. OISST V2 covers more than 30 years and provides a good chance to analyze SST variations at interannual and inter-decadal scales. It has been widely used in MHW-related scientific researches [25,26]. It is important to recognize that in spite of the fact that the majority of research, including the most recent global assessments, have characterized MHWs using OISST, this product may differ from in-situ observations such as Argo in certain regions, hence undermining its use for MHW analyses. For example, in tropical regions, with strong variability of heat fluxes and dense cloud cover, small differences can be observed between satellite-based and buoys, ships, or drifters’ data [27]. In the SCS, the lowest area-averaged SST value (26.4 °C) occurs in winter (December-February), a 1.7 °C drop from the annual mean of 28.1 °C (Figure 1c). The area-averaged SST value begins to rise in March, and it surpasses 29 °C from May to September, which is 1.0 °C (i.e., nearly 1 standard deviation) higher than the annual mean value. Thus, the five consecutive hottest months are considered the warm season in this region (Figure 1c). In contrast to other months, MHWs in these five consecutive hottest months indicate extraordinarily high SST conditions.

The impact of atmospheric conditions and related atmospheric circulation on MHWs trends were examined by composing the anomalies of associated atmospheric variables in the National Centers for Environmental Prediction/National Center for Atmospheric Research (NCEP/NCAR) reanalysis (http://www.cdc.noaa.gov (accessed on 27 September 2022)) [28]. These variables include monthly mean geopotential height (GHT), zonal and meridional 10 m wind, cloud cover, surface net downward solar shortwave radiation, and vertical velocity fields from 1948 to 2020. These data have a horizontal resolution of 2.5° × 2.5° and 17 vertical levels ranging from 1000 to 10 hPa. Note that the atmospheric circulation changes based on the ERA5 and NCEP/NCAR datasets have been compared in several previous articles, and we find that the results from two different reanalysis datasets are highly consistent [29,30]. We thus present the results based on the NCEP/NCAR dataset, which has longer coverage than ERA5.

Daily mean geopotential height at 500 hPa level (GHT500) from NCEP/NCAR Reanalysis 1 dataset was used to devise the composite maps (https://psl.noaa.gov/data/gridded/data.ncep.reanalysis.html (accessed on 27 September 2022)) [28]. The spatial coverage of data is 0.0°E to 357.5°E, 90.0°N to 90.0°S, with a horizontal resolution of 2.5° × 2.5°. The data ranges from 1948 to 2020, and daily climatology is derived from the 30 years (1983–2012) and smoothed with the method of a centered-moving mean of 31 days.

The National Climate Center (NCC) of the China Meteorological Administration (CMA) recommended four operational indices based on the 5880-gpm contour at 500 hPa to objectively describe the status of the Western North Pacific Subtropical High (WNPSH) in monitoring service, including area index (AI), intensity index (INI), ridgeline position index (RI), and Western ridge point index (WRI). The four WNPSH indices represent the WNPSH’s area, intensity, northern extension, and western borders, respectively. More details about the definitions of these indices can be found in Liu et al. (2012) [31]. In our work, the monthly WNPSH indices from 1982 to 2020 were obtained from the website of NCC (http://cmdp.ncc-cma.net/cn/index.htm, accessed on 10 April 2021).

Our analysis was limited to 1982–2020, and the baseline in our work for the climatology is 1983–2012, which is shifted by two years from the standard normal period of 1981–2010 because the satellite SST series began in 1981.

2.3. Identification of Marine Heatwaves (MHWs) and Indices

In most cases, two types of definitions based on absolute and relative thresholds are employed to identify the heat waves. For better comparisons with other regions and the globe, we adopted a qualitative definition proposed by Hobday et al. (2016) [32]. According to the definition, an MHW is an unusually warm discrete event lasting for five or more days and throughout which the daily mean SST exceeds a relative threshold. Since MHWs are defined as the occurrence of daily SST anomalies at the upper tail of the probability distribution, a higher threshold, such as a 90th percentile, is employed, which can detect a relatively large number of MHWs during the warm season. Here, the 90th percentile of a 31-day moving-centered window from 1983 to 2012 establishes the relative threshold for each calendar day during the warm season (MJJAS). Two consecutive MHWs with a two-day or less gap were treated as a single event. Once the MHW events at any grid point are detected, a set of summary statistics can be derived for MHWs, including, for example, frequency, duration, and intensity. For an MHW, duration (D) is the period between the start (ts) and end (te) date of the MHW event. Considering the cumulative effect of thermal stress on marine ecosystems, we utilized the cumulative intensity to represent the MHW intensity. Cumulative intensity is the integral of intensity over the event duration (units: °C days). It is equivalent to previously used metrics such as Degree Heating Week, a method utilized by the NOAA Coral Reef Watch program for real-time coral bleaching alerts [33]. In the following analysis, for a particular grid point in the warm season of a given year, MHWN represents the number of MHWs in this year; MHWD represents the average duration of all MHWs in this year, and MHWI represents the average cumulative intensity of all MHWs in this year. All of the MHW indices listed above are detailed in

Table 1. Definitions of annual MHW indices. $N$ is the number MHWs in the warm season of each year; D is the period between the start and end date of an MHW event; ${\int }_j^D(SST(t_j)-SS{T}_m(t_j))dt$ is the cumulative intensity which is calculated as the sum of the difference between the daily SST for a certain day and the daily climatological mean ($SS{T}_m$) over the MHW duration (D). Notably, here the anomaly refers to a deviation of daily SST from the climatological mean which we used for the fixed-baseline period (1983–2012).

Index	Definition	Symbol or Formula	Unit
MHWN	Number of MHWs in the warm season in each year	N	Count
MHWD	Average duration of all MHWs in the warm season of each year	${\displaystyle \sum _{i=1}^N}\left(D_i\right)/N$	Days
MHWI	Average cumulative intensity of all MHWs in the warm season of each year	$\left({\displaystyle \sum _i^N{\displaystyle {\int }_j^{D_i}(SST(t_{ij})-SS{T}_m(t_{ij}))}dt}\right)/N$	°C•days

. The code of the index calculation had been made available as a free heatwaveR package [34].

2.4. Regional Marine Heatwaves (MHWs) and The Relevant Circulation Anomalies

A regional MHW was determined when the regional averaged (5°S~25°N and 100°~125°E) daily mean SST was above the 90th percentile value and persisted for five or more days. Here, 37 regional MHWs were identified in the SCS from 1982 to 2020. WNPSH is a large-scale and persistent high-pressure circulation system that would last for several weeks. To avoid the influence of the previous MHW on subsequent ones, we selected 28 regional MHWs that were more than 7 days apart from the end day of the previous MHWs. Then, the averages of the 500-hPa geopotential height of the 28 MHWs were counted by giving each event an identical weight to find the relevant circulation anomalies.

2.5. Empirical Orthogonal Function (EOF) Analysis

We used an empirical orthogonal function (EOF) analysis to investigate the spatiotemporal patterns of MHWN, MHWI, and MHWD of these detected MHWs in order to emphasize the key modes responsible for the interannual fluctuation in MHW.

The EOF method aims at decomposing a data set into a product of a set of spatial patterns [35] and time series as follows:

$$d\left(x_i,t_j\right)=\sum _{k=1}^M{\varnothing }_k\left(x_i\right)a_k\left(t_j\right)$$

(1)

where $d$ is the data set, $x_i$ is the position, $t_j$ is the time, ${\varnothing }_k$ is the spatial pattern of the $k-th$ mode, $a_k$ is the amplitude time series of the $k-th$ mode, and $M$ is the number of spatial elements. The functions ${\varnothing }_k$, which are the basis functions, are chosen to be orthogonal to the other basis functions because the information would be repeated if not orthogonal, and the unit module is a common choice with no loss of generality:

$$\sum _{k=1}^M{\varnothing }_i\left(x_k\right){\varnothing }_j\left(x_k\right)={\delta }_{ij}$$

(2)

where ${\delta }_{ij}$ is the Kronecker delta.

Each of the EOF patterns identified has a corresponding time series showing how the amplitude of the pattern changes in time, known as the principal component (PC). The PC time series may be obtained by projecting the EOF pattern onto the original field at each time step to determine the sign and amplitude of the pattern at any given time. A given mode can be reconstructed by multiplying the EOF (space) by its PC (time).

The significance of the EOF modes was determined using North statistical testing [36]. That is, computing their sampling error as $\lambda (2/N{)}^{1/2}$, where $\lambda$ is a given eigenvalue and N the number of realizations, here N = 1200.

2.6. Testing for Significance of Correlations and Trends

When testing the significance of the correlation between two time series or the presence of trends in a time series, one important question to consider is how large sample correlations and sample trends can be, even if the stochastic processes, which generate the series, are not correlated and do not exhibit any trends. Firstly, we need to assume that processes X and Y share no correlation or that segments of length L of the process have no trend. Standard procedures are available in the literature, namely p-value for correlations and Mann–Kendall for trends [37], that there are “no correlations” between the underlying processes, and trends can hardly appear in limited segments of an infinite stationary time series.

In the case of correlations, the assumption is that the underlying processes are stationary (free of systematic trends) and serially independent, i.e., $X_t$ and $X_{t+1}$ for any t are independent. In the case of trends, the assumption is the independence of $X_t^{\prime }$. However, in geophysical cases, these assumptions are not satisfied. The result is that the null hypotheses are often falsely rejected (i.e., in cases where there are no correlations or no trends) than stipulated by the significance level (normally 5%).

2.7. Mann–Kendall Abrupt Change Test

Abrupt change occurs when the climate system is forced to cross a threshold and denotes a fast transition from one state to another [38]. The Mann–Kendall abrupt change test, also known as the sequential Mann–Kendall–Sneyers test, is one of the most effective methods for determining when an abrupt change occurred in a time series [39]. It has been widely applied to change point detection for long-term time series data (e.g., hydrological and meteorological climate changes). Our study employed this approach to examine the change points in MHWs.

For the time series $x_n$ (n is the length of the data set), the order series ($S_k$) is given as follows:

$$S_k={\displaystyle \sum _{i=1}^k{r}_i,(k=2,3,4,\dots ,n)}$$

In which,

$$r_i=\left\{\begin{array}{l}+1,x_i>x_j\\ 0,x_i\le x_j\end{array}\right.,(j=1,2,3,\dots ,i)$$

where $x_i$ and $x_j$ are the sequential data values, and the statistic ($U{F}_k$) is defined as:

$$U{F}_k=\frac{\left[S_k-E(S_k)\right]}{\sqrt{Var(S_{{}_k})}},(k=1,2,3,\dots ,n)$$

where $U{F}_1=0,$ $E(S_k)$ and $Var(S_k)$ are the average value and variance of $S_k$, which can be calculated by the following equations:

$$E(S_k)=\frac{n(n+1)}{4}$$

$$Var(S_k)=\frac{n(n-1)(2n+5)}{72}$$

Then, the $U{B}_k$ is calculated by repeating the above process in the order $x_n,x_{n-1},\dots ,x_3,x_{2,}{x}_1$ of the time series which makes $U{B}_k=U{F}_k(U{B}_1=0,k=n,n-1,\dots ,3,2,1)$. In this study, we set p values as 0.05. The intersection points of the $UB$ and $UF$ curves indicate the possible turning year of the trend within a time series [40]. The trend turning point is considered significant at the corresponding threshold level (i.e., ±1.96 for a 95% significance level).

3. Results

3.1. Climatology and Long-Term Trends in Annual MHWs Indices

Figure 2 displays the spatial distributions of the multiyear warm season average MHW indices in the SCS. Notably, MHWN, MHWD, and MHWI are not spatially uniform in the region. More common, longer, and more intense MHWs primarily occur in northern SCS with shallow depth, particularly around the Beibu Gulf and the offshore areas of the Pearl River estuary (Figure 2a–c). Longer-lasting MHWs also occur at a core section of the deep basin, where the average depth is over 4000 m. MHWs are less frequent in southern SCS, despite the fact that the climatological mean and 90th percentile SSTs of MJJAS in this region are both exceeding 30 °C, which is higher than in other locations (Figure 2a–c).

During the study period (1982–2020), the MHWs in the SCS was enhanced in terms of increasing frequency, prolonged duration, and strengthening intensity (Figure 2d–f). The strongest signals of the enhancements are primarily found in coastal regions at shallow water and higher latitudes in the north SCS and the enhancements are almost twice those in the south SCS. In addition, though the linear trends of MHWN in the central SCS were insignificant, the trends of MHWD and MHWI were both significantly positive (Figure 2e,f). It means that MHWs in the central SCS have been getting even substantially severer (i.e., longer duration and/or high intensity). Thus, with the largest number of coral reefs, the Nansha Islands located in the central SCS face the fiercest threat from MHWs.

3.2. Dominant Modes of Marine Heatwave Variations

To investigate the spatial structures and temporal properties of the MHWN, MHWD, and MHWI in the SCS, EOF analysis is performed (Figure 3, Figure 4 and Figure 5). The first leading EOF (EOF1) modes of MHWN, MHWD, and MHWI explain 56.8%, 64.5%, and 52.6% of the total variance, respectively. The EOF1 modes are well separated from others according to the North test. All of the EOF1 spatial patterns have a homogeneous structure over the SCS, and the corresponding principal components (PC1s) time series show clear increasing trends (Figure 3, Figure 4 and Figure 5). In addition, there are notable interdecadal variations in the MHWN, MHWD, and MHWI time series, with a turning point occurring around the year 2010 based on the Mann–Kendall test method and the rapid rises since then (Figure 6).

The remaining spatial patterns describe a small percentage of the variance (i.e., approaching or below 10%), dominated mainly by local features. Second EOF (EOF2) modes of MHWN, MHWD, and MHWI explain 8.1%, 8.4%, and 8.1% of the total variations, respectively. The spatial distribution variations of these modes are the meridional dipole pattern with anomalies of opposing signs in the north and south of the SCS (Figure 3, Figure 4 and Figure 5). The corresponding PC2 time series exhibited high interannual variations. It is worth noticing that the amplitude of PC2 has increased dramatically since 2010, implying that the dipole pattern of frequency, duration, and intensity of the MHWs has been considerably enhanced, which is consistent with the earlier findings [41]. The spatial patterns for the third EOF (EOF3) modes, which explain only 5~6% of total variations, exhibited a triple pattern with notable positive anomalies centered at the center SCS. PC3 time series exhibited an opposite trend to PC1 and PC2, and the triple pattern has become less obvious since 2010.

Changes in extremes can be linked to changes in probability distribution’s mean, variance, shape, or all of these [42,43]. As mentioned above, MHWs have noticeable decadal variations in frequency, duration, and intensity, with more frequent, longer, and severer MHWs in the recent decade. Figure 7 depicts a statistical analysis of the regional averaged daily SSTAs in the warm season of the four decades. It is obvious that area-averaged SSTAs have greatly increased from the first decade, 1982–1990 (−0.24 °C), to the last decade, 2011–2020 (0.25 °C), with a significant warm toward-shift, up to 0.49 °C (Figure 7). In addition, the histogram distribution has been substantially wider and flatter in the recent decade (kurtosis = 0.77 in 2011–2020, kurtosis = −0.27 in 1982–1990). The shape in the last decade has become skewed towards the hotter part (skewness = 0.18 in 2011–2020, skewness = 0.10 in 1982–1990), and the dispersion degree is greater. Increased occurrences of MHWs, with stronger intensity and longer duration, are a result of such a change in distributions. According to the findings, the mean temperature and natural variability both contribute to the incidence of MHWs in the SCS. In other words, variations in MHWs cannot be explained solely by mean SST warming from anthropogenic forcing [44].

3.3. Connection with Large-Scale Atmospheric Circulation

3.3.1. Abnormities of Atmospheric Circulation

MHWs can be modulated by various mechanisms, including localized air-sea heat flux, atmospheric and oceanic conditions, and remote large-scale teleconnections, according to previous research [4,45]. Given that atmospheric conditions were the primary driver of MHW priming [46,47], we investigate the impact of atmospheric conditions and related circulations on MHWs trends in this section. For 1982–1991 (an early period with low-frequency MHWs) and 2011–2020 (a recent period with high-frequency MHWs), composite analysis of the atmospheric and near-surface meteorological fields, namely 10 m winds, total cloud, surface net down-ward solar shortwave heat flux and vertical velocity anomalies, are calculated, respectively.

The distributions of atmospheric and near-surface field anomalies in 1982–1991 and 2011–2020 are shown in Figure 8, respectively. From 2011 to 2020, a notable anomalous anticyclone circulation was observed in the SCS, with the center at 116°E, 13°N, while anomalous cyclonic circulation existed in the SCS from 1982 to 1991. The anticyclone circulation in the recent decade has the potential to create an abnormal easterly and northeasterly flow from the center to the south, including the southern coast of Vietnam. The unusually strong easterly wind may diminish offshore Ekman transport and the upwelling in the east of Vietnam, weakening the upwelling cold core (110°E, 12°N) [23], increasing the surface temperature of central SCS, as compared to the early decade. Except for the southern part of SCS, the SCS is accompanied by decreased cloud cover and increased solar radiation reaching the surface, which is likely to be the main reason for the more frequent MHWs from 2011 to 2020 (Figure 8d,f). The longitude-altitude cross sections of vertical velocity anomalies along 12.5°N for 1982–1990 and 2011–2020 are shown in Figure 8g,h. From 2011 to 2020, weak sinking motion dominated the entire air column above the center of the SCS, with the largest sinking motion occurring between 500 and 300 hPa, while rising motion occurred over the SCS for the early period, with the largest rising motion appearing at the lower troposphere (Figure 8). In the recent decade, the anomalous descending air motion could have been attributed to the anomalous anticyclone, which led to less cloud, and more solar radiation reaching the surface.

3.3.2. Western North Pacific Subtropical High (WNPSH)

As detailed in the preceding subsection, these anomalies, including decreases in wind speed, more solar radiation, less cloud cover, and stronger downward air motion, are likely linked to the WNPSH’s westward shift and strengthening. WNPSH is the western extension of the North Pacific Subtropical High, which extends to the East Asian coast. Stronger WNPSH generates more floods and heat waves in East Asia on interannual to interdecadal timescales [48,49].

Composite patterns of the GHT500 anomalies and the mean position of WNPSH, which is denoted by the 5870 gpm contour line at 500-hPa, for the first and last decade are shown in Figure 9, respectively. For the recent decade (2011–2020), notable positive anomalies appeared throughout the entire SCS, while the opposing distributions of the GHT500 anomalies were in the first decade (1982–1990) (Figure 9). Moreover, Figure 9b shows a significant westward extension of the WNPSH compared with the climatological position (green and blue dashed lines). These findings suggest that the intensification and westward extension of the WNPSH may be important in establishing MHWs in the SCS, as they lead to strong subsidence, decreased cloud cover, and more solar radiation.

It is worth mentioning that, since the 1980s, the SCS has become warmer, accompanying global warming [50] (Figure S1, in Supporting Information S1). Previous research has also indicated that long-term warming caused by human activities has strengthened the WNPSH [51,52]. Before examining the MHW changes connected to WNPSH, linear trends of the PC1 time series and GHT500 should be removed (i.e., removing the linear trend from the time series and remaining natural variability). Regression analysis is a powerful statistical method that examines the relationship between two or more variables of interest. Figure 10 displays the patterns of detrended GHT500 regressed on the detrended PC1 time series of MHWN, MHWD, and MHWI. The whole SCS was dominated by positive GHT500 anomalies for all the detrended PC1s time series, with the positive maximum over the north SCS (15–25°N, 105–125°E, black rectangle in Figure 10). The observations support prior findings that an unusual high-pressure system was conducive to more heat waves in Southeast Asia [49]. The strengthening and westward expansion of the WNPSH could result in strong descending motion and solar radiation, both of which may have contributed significantly to the region’s sea surface warming.

Correlation coefficients between MHW and WNPSH indices from 1982 to 2020 were calculated before and after detrending (

Table 2. Correlation coefficients (r) between the WNPSH indices and the MHW indices temporal series from 1982 to 2020, with and without detrending.

Index	r			r (Detrended)
	MHWN	MHWD	MHWI	MHWN	MHWD	MHWI
Area index (AI)	0.76 *	0.76 *	0.75 *	0.70 *	0.67 *	0.66 *
Intensity index (INI)	0.7 5 *	0.77 *	0.75 *	0.69 *	0.69 *	0.68 *
Ridge line index (RI)	−0.27	−0.19	−0.20	−0.36	−0.28	−0.27
Westward ridge point index (WRI)	−0.75 *	−0.64 *	−0.63 *	−0.57 *	−0.51 *	−0.51 *

Asterisk * denotes that coefficients are significant and exceed the 0.01 level (0.39) under the t-test [53].

). The absolute correlation coefficients (r) between the area index, intensity index, and Westward ridge point index of the WNPSH and PC1s time series before and after detrending were all greater than 0.50 (exceeding the 0.01 confidence level), although the correlation coefficients decreased after detrending. However, there is no clear link between the RI, representing the meridional displacement of WNPSH and PC1s, with r below 0.05 confidence level. Both the long-term trend and the variabilities of the zonal displacement and intensity of WNPSH are likely responsible for the MHWs, according to the above findings. Regionally, the results were congruent with those of other studies on land heatwaves [49].

To get a deeper comprehension of the role of the WNSPH in the genesis of MHWs, we further examine the evolution of the synoptic characteristics connected to MHWs. The composite charts from 6 days before the starting date of 28 MHW events selected from all 37 MHWs (see Subection 2.4) were placed at a 2-day interval from top to bottom (Figure 11). These patterns illustrate the temporal evolution of GHT 500 anomalies linked with MHWs in the SCS. The positive anomaly center of GHT 500 emerges around (135°E, 15°N) 6 days before the MHW onset (Day-6) (Figure 11a). This positive center steadily strengthened, traveled westward over time, and was positioned at the center of SCS (122°E, 15°N) on 2 and 4 days before the MHWs onset (Day-2 and Day-4, Figure 11b,c). On Day 0, when the MHWs start, these anomalies greatly intensify, spread southwestward, and finally cover most of the SCS. According to the findings, there is considerable correspondence between the temporal evolution of the GHT500 anomalies and MHW emergence. It confirmed that the westward displacement and intensification of the WNPSH may be crucial for the formation of MHWs in the SCS. Previous studies of the WNPSH motion also indicated that the WNPSH mechanism played a vital role in the onset and termination of terrestrial heat waves in southern China [49].

4. Discussion

Human-induced climate change is causing upper ocean warming globally and has led to observed increases in the frequency, length, and severity of MHWs in recent decades [5,11]. These MHWs threaten marine biodiversity and ecosystems, increase the likelihood of extreme weather, and negatively impact the fisheries, aquaculture, and tourism industries and therefore attract significant scientific and public interest [4,54]. Little work has been conducted on SCS MHWs, despite the frequency, length, and severity being comparable to other typical areas, such as the Mediterranean [6] and the Tasman Sea [12].

It is essential to acknowledge that there are still gaps in our understanding of MHWs. From the definition of MHWs [32], it can be seen that the baseline period of SST is crucial and can directly affect the percentile metrics of MHWs. Hobday et al. suggested a fixed-baseline period of 30 years, considering the limited time coverage of SST data and the influence of natural oscillation such as ENSO. The method of defining MHWs using the relative threshold of fixed-baseline periods has been broadly adopted. However, as the ocean warms, the real temperature itself will shift, and rising mean SST might increase SCS MHWs. Oliver (2019) [43] emphasized that SST mean warming, not variability, dominated the MHW trend. Therefore, some researchers propose to use a moving baseline to describe MHWs, which is more advantageous for describing short-term SST fluctuations [4,55]. Both methods of calculating baseline climatologies are complementary; the choice depends on the questions posed. A fixed-baseline period is particularly relevant when examining ecosystems with a slow or even no adaptive capacity and displaying overall changes, including long-term warming effects. Moving-baseline periods are most beneficial for considering ecosystems with a faster adaptation capacity or greater mobility and for studying physical climate characteristics with a focus on variability rather than long-term warming, as this approach primarily reflects shorter-term variability changes. To this end, research suggested that combining both approaches can effectively distinguish the effects of long-term warming from changes in the magnitude of variability [4].

The primary focus of this work is on the first dominant modes of MHWs, and their connections to the major atmospheric circulation system—WNPSH. Multiple processes from ocean basins affect the variability of the WNPSH [56,57]. Previous studies have demonstrated that the unusual WNPSH in summer is closely related to the ENSO decaying phase, with the rapid decay of El Niño being more favorable to a stronger WNPSH than the slow decaying of El Niño. If a La Niña event has already formed in summer during the El Niño decaying year, the strongest WNPSH will occur [13,58]. An incipient La Niña suppresses central Pacific local convection, generating a westward propagated atmospheric Rossby wave that could strengthen the WNPSH [56,57]. Enhanced convection over the Maritime Continent associated with the La Niña strengthens the WNPSH by modulating local Hadley circulation with increased anticyclonic circulation over the northwestern Pacific [59].

In addition, it is necessary to emphasize that MHWs are modulated by a wide range of physical mechanisms, not just atmospheric circulation, spanning from local dynamics to distant processes via teleconnections on various spatial and temporal scales [42,45]. In recent studies, researchers highlighted the significant interannual and interdecadal variabilities of SCS MHWs events and their associations with large-scale climate modes, such as ENSO [14,60]. They note that more MHWs were likely to occur in the SCS during the following summer after the El Niño events, which has important implications for the predictability of the SCS MHWs. Moreover, some research suggested that other factors, such as local upwelling and entrainment, which were weakened or even disappeared, and the Arctic amplification links to the Rossby wave, which moved more slowly in the upper-level flow, could enhance the possibility of MHWs [13,61]. Yet, in comparison with the AHWs occurring over the coastal hinterland of the South China [49], there is far less knowledge on which factors and to what extent could cause the interdecadal and interannual variability in MHWs in the SCS. A deeper understanding of the mechanisms of MHWs and the relative effects of these drivers in the SCS requires more investigation, especially some quantitative studies.

Recent MHWs with high scientific and public attention include the 2014/15 North Pacific MHW, the 2015/16 and 2017/18 Tasman Sea MHWs, etc. [12,25,46]. These events heavily impacted fisheries and aquaculture, in addition to causing coral bleaching and seagrass mortality. Much less concern is given to the severe MHWs in the SCS, let alone some mild or moderate events, despite their devastation of marine ecosystems and the services they provide. On the Nansha Islands in the SCS, there are around 8000 km² of coral reefs that provide a substantial number of people with goods and services. The multiple metrics (duration, intensity and heating rates, etc.) and thermal structure of MHWs that induce coral bleaching in this region have not been previously studied, though there are notable relationships between MHWs and coral bleaching. Finding and exploring the connection between cumulative intensity (i.e., accumulated heat stress) and coral mortality at regional scales is highly encouraged. Moreover, the local government has not yet implemented any systematic real-time monitoring, sub-seasonal, and seasonal forecasting systems to help marine managers and industries cope better with MHWs. Consequently, specific studies will be required to solve these issues in the SCS.

5. Conclusions

MHWs, extreme water events, can occur at any time of the year. This paper mainly focused on the warm season MHWs in the SCS, which can easily exceed organisms’ thermal thresholds (e.g., coral reefs) and lead to severe coral bleaching and mass mortality. Results showed that during the study period (1982–2020), MHWs had been enhanced across most of the SCS in terms of increasing frequency, prolonged duration, and strengthening intensity. The strongest signals of the enhancements were primarily found in coastal regions at shallow water and higher latitudes in the north SCS, and the enhancements were almost twice those in the south SCS.

The multiscale variations of warm season MHWs metrics (i.e., MHWN, MHWD, and MHWI) in the SCS were analyzed using an EOF method. The first leading EOF (EOF1) modes of MHWN, MHWD, and MHWI explain 56.8%, 64.5%, and 52.6% of the total variance, which account for more than half of the total variance. All of the EOF1 spatial patterns reveal a consistent variation of a homogeneous structure. The corresponding PC1s time series display distinct interdecadal variations with a turning point in each time series around 2010. Second EOF (EOF2) modes of MHWN, MHWD, and MHWI explain only 8.1%, 8.4%, and 8.1% of total variations, respectively. The spatial distribution variations of these modes are the meridional di-pole pattern with anomalies of opposing signs in the north and south of the SCS. The corresponding PC2 time series exhibited high interannual variations. It is worth noticing that the amplitude of PC2 has also increased dramatically since 2010.

Understanding the variability of MHWs and their relation to large-scale circulation can be important for predicting MHWs. By linking the first dominant modes with the large-scale atmospheric circulation, it is clearly shown that the warm season MHWs in the SCS were commonly featured with anomalous high-pressure system-WNPSH, and accompanied by local calm wind, descending motion, less cloud cover, and enhanced solar radiation heating at the surface. In the majority of the SCS, the regression coefficients between the natural variabilities of the three PC1 time series and the natural variabilities of GHT 500 are all significantly positive. Our study highlights the more important role of WNPSH’s area, intensity, and westward ridge point location than ridge line location. The evolution of the MHWs was associated with the westward displacement of WNPSH, with the notable and amplified positive anomaly in GHT 500 migrating westward.

MHWs have become one of the major challenges facing the sustainable development of marine ecosystems and environment. In the context of continuous global warming, MHWs are predicted to increase in frequency, duration, and intensity in most regions [3,5,6]. Regional capacity building for MHW monitoring, prediction, and forecasting must be strengthened. Recent research demonstrated that the occurrence of MHWs was mostly related to large-scale atmospheric and oceanic modes, especially low-frequency climate modes (such as ENSO and WNPSH). Therefore, we highlight the potential of predicting and forecasting and provide a foundation for a much-needed operational MHW forecast system [55]. There is still work to be done in the future to improve the accuracy of the prediction system for MHWs through more research on the synergies among various climate mode impacts on MHWs, the development of high-resolution models, and artificial intelligence prediction technologies over the SCS.