The Trend and Interannual Variability of Marine Heatwaves over the Bay of Bengal

: Marine heatwaves (MHWs) are long-lasting extreme oceanic warming events that can cause devastating effects on warm-water corals and associated ecosystems. The linear trend and interannual variability of MHWs over the Bay of Bengal (BOB) during 1982–2020 are investigated by a high-resolution daily sea surface temperature (SST) dataset. In regions where warm-water coral reefs are concentrated, annual MHW days and frequency signiﬁcantly increase during 1982–2020, at rates exceeding that of the global mean. The coldest boreal winter season witnesses signiﬁcant and steady increase trends in MHW days and frequency. In contrast, the trend is insigniﬁcant in the climatological warmest season (March to June) south of 15 ◦ N in the BOB, mainly due to large interannual variability. El Niño and Southern Oscillation (ENSO) dominates the interannual variability of BOB MHWs, which are highly consistent with the evolution of the mean SST. The negative phase of North Atlantic Oscillation (NAO) also modulates the occurrences of MHWs, especially over the northeastern BOB. The two climate modes synergistically explain about 50~70% of the interannual variances in the BOB’s MHWs. Correlation analysis reveals that south of 15 ◦ N in the BOB, the effect of El Niño on MHWs is evident from the boreal autumn of its developing phase to the boreal summer of its decaying phase, along with limited inﬂuence from NAO. However, in the northeast of the BOB, the effect of El Niño merely emerges from April to August of its decaying stage. In comparison, boreal winter-to-spring NAO exerts a strong control over March-to-June MHWs in the northeastern BOB. The results suggest that various climate modes may jointly or separately inﬂuence MHWs at certain seasons and locations, which is important for the seasonal prediction of MHWs. Indeed, when combining the Niño3.4 mature winter index and boreal winter-to-spring NAO index to build a regression model, it is more effective in reproducing the BOB’s MHW frequency compared to the Niño3.4 index alone. despite the exclusion of ENSO’s effect. The results are similar if we replace the


Introduction
Marine heatwaves (MHWs) are extreme sea surface temperature (SST) warming events lasting for several days, weeks or even months. The resultant coral bleaching, seagrass meadow destruction and species migration are devastating for marine ecosystems and fishery production [1][2][3][4]. However, MHWs receive great attention only after the 2011 West Australian MHW, despite their devastating effects [5,6]. Under global warming, MHWs are projected to occur more frequently with prolonged days [7][8][9], further damaging or even destroying the living conditions of marine life. This highlights the importance of investigating the drivers and mechanisms of MHWs, and hence improving the short-term predictions and long-term projections [10][11][12].
Previous studies mostly focus on the mechanisms of MHWs in a few hot regions or discuss the mechanisms on a global scale. The mechanisms of individual MHW cases display large differences across regions, but mainly relate to local oceanic and atmospheric conditions, air-sea coupling and climate mode-associated remote teleconnections [6,[13][14][15]. For example, the anomalous air-sea heat flux generated by the atmospheric high-pressure system lead to the 2003 Mediterranean MHW [16], and abnormal ocean advection causes the 2015/16 Tasman Sea MHW [17] and the southeastern tropical Indian Ocean MHW [18]. Climate modes play an important role in modulating the variability of MHWs, but the dominant modes substantially vary across regions. Generally, MHWs are closely associated with the El Niño and Southern Oscillation (ENSO), the Indian Ocean basin mode (IOB) and dipole mode (IOD) and the Pacific Decadal Oscillation (PDO) in Indo-Pacific regions [19][20][21][22][23][24][25][26]. The North Atlantic Oscillation (NAO) is also shown to be connected with the formation of Indian Ocean MHWs [13,24].
The Bay of Bengal (BOB) is a critical region of the Indian Ocean because it supports large marine primary productivity and is strongly coupled with the Indian monsoon [27]. Additionally, the BOB suffers powerful and deadly tropical cyclones as it supplies warm SSTs and high humidity year-around. Additional extreme and prolonged warming beyond the high SSTs can cause severe disasters in marine ecosystems, such as the 2010 Andaman sea MHW [28], and may even regulate the Indian summer monsoon [24] and cyclone activities. However, the BOB MHWs during the last decades have lacked systematic and in-depth investigations, with their evolution features in recent decades remaining unclear. For example, El Niño is revealed to have caused two peak warming in the North Indian Ocean during its developing phase and decaying phase, respectively [29][30][31]. Therefore, the BOB MHWs may also be more active during the North Indian Ocean peak warming periods. The detailed relationship between the BOB and climate modes has not yet been systematically revealed. Global warming is supposed to cause an increasing trend in the duration and frequency of MHWs over most oceans [13,15]. Given that the warming trend is predicted to continue in the near future, background condition changes may further fuel the BOB MHWs when favorable climate modes occur.
Therefore, the present study investigates the linear trend and interannual variability of BOB MHWs during 1982-2020, when satellite observation data was widely available. We show that the BOB's annual MHW days and frequency significantly increase at rates that are larger than that of the global mean. Causality analysis shows that the increasing trends of BOB MHWs are clearly driven by the rise in atmospheric CO 2 concentration. However, the linear trend is not significant in the warmest season, especially for boreal spring. The interannual variability of BOB MHWs is large and mainly dominated by ENSO, with the negative phase of the boreal winter-to-spring NAO also being effective in regulating their occurrences, especially over the northeastern BOB. The two climate modes synergistically explain about 50~70% of the interannual variances in the BOB MHWs.
The rest of this paper is organized as follows: Section 2 describes the definition of MHWs, data and methods; Section 3 shows the climatology and linear trend of the BOB MHWs' properties; Section 4 investigates the interannual variability of the BOB MHWs and underlying processes; and Section 5 presents the discussion and summary.

Data
Coral reef data (WCMC008 v4.1) shown in Figure 1a is collected from 1954-2018 by the UNEP World Conservation Monitoring Centre (UNEP-WCMC) and the WorldFish Centre, in collaboration with the WRI (World Resources Institute) and TNC (The Nature Conservancy). Data sources include the Millennium Coral Reef Mapping Project [32,33] and the World Atlas of Coral Reefs [34]. It shows coral reef distributions in the tropical and subtropical regions and is the most comprehensive empirical global dataset of warmwater coral reefs to date. The high-resolution daily SST analysis (Optimum Interpolation SST, OISST) are developed using the optimum interpolation method [35]. The analyses have a spatial grid resolution of 0.25 • and a temporal resolution of 1 day, and the cold bias against Argo is about −0.14 • C on a global average and −0.28 • C in the Indian Ocean during 2016-2019. It is the most used dataset for the detection of MHWs. We also utilized it to detect MHWs from January 1982 to December 2020 in the BOB. Atmospheric carbon dioxide measured at Mauna Loa provides the longest record of atmospheric CO 2 concentrations [36,37]. The monthly mean NAO index, since 1950, is provided by the National Weather Service Climate Prediction Center.

MHW Definition
Following Hobday et al. (2016), an MHW is defined as a prolonged discrete anomalously warm-water event in which SST warming exceeds a threshold of the 90th percentile based on the climatological daily mean and lasts for at least 5 days. The 90th percentile and climatology are calculated within an 11-day window centered on each calendar day, and then smoothed in a 31-day running mean window by the daily SST from 1982 to 2020 [14,38,39]. By this definition, MHWs can be identified at any time of the year, including the local winter season. An MHW is primarily measured by the metrics of its days (MHW days, number of days between the specific period), frequency (MHW frequency, number of MHWs during a specific period) and mean intensity (MHW intensity, averaged SST anomaly during each MHW event). In this study, MHW days and frequency are attributed to be zero when no MHW occurs during a certain period.

Causality Analysis
Following Liang (2014) [40], the causality between two time series (X 1 and X 2 ) can be judged by the Liang-Kleeman information flow theory. The exchange of information between two events indicates the quantity as well as the direction of the causality. Liang (2014Liang ( , 2015Liang ( , 2016 [40][41][42] show that if the development of X 1 is independent of X 2 , then there is zero information flow from X 2 to X 1 . For the linear systems, the maximum likelihood estimator of it can be estimated as: where C ij is the sample covariance between X i and X j , C i,dj is the sample covariance between X i and X j , which is defined as X j = { X j,n+1 − X j,n /∆t}, and ∆t is the time step. If |T 2→1 | = 0, X 2 does not cause X 1 ; if |T 2→1 | = 0 and the value passes the statistical significance test, X 2 is deemed to cause X 1 [40][41][42]. Stips et al. (2016) suggest that if |T 2→1 | |T 1→2 |, X 2 is also supposed to cause X 1 ; if they are comparable and both pass the significance test, there is a reciprocal causation. The calculation script of the method is available at http://www.ncoads.cn/article/show/63.aspx, accessed on 8 March 2022. This method has been widely applied in different studies to explore the causality between different time series [43][44][45][46][47].

Study Area
The devastating nature of extreme warming on coral polyps results in warm-water coral reefs being increasingly threatened by MHWs. Warm-water coral reefs are highly restricted in areas of warm, shallow and clear water to produce the copious quantities of limestone necessary for reef formation. To focus on the influence of MHWs on coral bleaching, we selected three regions in which coral reefs are concentrated, that is Region 1 (77 • -83 • E, 5 • N-12.   than the other two regions. Ocean temperature varies seasonally following the migration of the subsolar point. In climatology, Regions 1 and 2 generally experience high SSTs above 28.5 • C from boreal spring to autumn (Figure 1e-g), with the highest SST appearing during boreal spring. The SST in Region 3 displays two distinct peaks in May and October, respectively, and is above 28.5 • C from April to November. Figure 1 displays the warm-water coral reef distribution and climatology of MHWs' annual days, frequency and mean intensity over the BOB for 1982-2020. In the Andaman Sea, where most coral reefs are located, MHWs take place more than 22 days per year in climatology ( Figure 1b). North of 15° N or east of 90° E, MHWs display relatively more days and exhibit a larger intensity than those in the western and southern regions of the BOB in climatology (Figure 1b,d). For the frequency (Figure 1c), MHWs occur around 2 times/year in most regions and display weak spatial variations. To illustrate the influence of MHWs on warm-water coral reefs, three regions were selected as the focus of this study (black box in Figure 1a), with Region 1 (77°-83° E, 5° N-12.5° N) in the southwest of the BOB, Region 2 (91°-100° E, 5° N-14° N) in the Andaman Sea and Region 3 (91°-95° E, 15.5° N-21° N) in the northeastern BOB. Region 1 displays much fewer MHW days (below 20 days per year) than the other two regions. Ocean temperature varies seasonally following the migration of the subsolar point. In climatology, Regions 1 and 2 generally experience high SSTs above 28.5 °C from boreal spring to autumn (Figure 1e-g), with the highest SST appearing during boreal spring. The SST in Region 3 displays two distinct peaks in May and October, respectively, and is above 28.5 °C from April to November.

Linear Trend of the BOB MHWs
We respectively calculated the area-weighted-mean MHW days, frequency and intensity of global oceans, Region 1, Region 2 and Region 3 (

Linear Trend of the BOB MHWs
We respectively calculated the area-weighted-mean MHW days, frequency and intensity of global oceans, Region 1, Region 2 and Region 3 ( MHW intensity that are not significant or negligible (Region 1 and Region 2), the values of the two-way information flow both fail to pass the significance test (Figure 2j,l). This also confirms the effectiveness of the information flow analysis in revealing the causality. Therefore, the significant increase trends in MHWs are clearly caused by the increase in atmospheric CO 2 concentrations (Figure 1a blue line). However, MHW intensity is above 1 • C in climatology and displays insignificant or even weak negative trends in three regions (Figure 2j-l), contrasting with the slightly increased global mean MHW intensity (Figure 2i). By definition, MHWs are extreme oceanic warming events with large intensity of SST anomalies. MHW intensity is calculated as the mean SST anomalies relative to the climatological daily mean of all the MHWs during specific time intervals, indicating that MHW intensity may also be high during early periods. As MHWs become more frequent and prolonged, the extended MHW days may decrease the averaged intensity of temperature anomalies for MHWs, causing negative trends in the MHW intensity in Region 2 despite the SST increase. about 1.45~1.6 times for the 1982-2000 mean in the three regions. The MHW frequency displays a trend that is weaker and insignificant (larger and significant) than the global oceans in Region 1, Region 2 and Region 3. Generally, MHWs display trends in days and frequencies that are larger in all three regions (exceeding 0.065 times/year) than global oceans (0.05 times/year). Causality analysis reveals that there is a significant information flow from the atmospheric CO2 concentration (Figure 1a blue line) to MHW days. In contrast, the information flow from MHW days to CO2 is negligible and insignificant. The results are similar for the information flow analysis between MHW frequency and CO2. Moreover, for the linear trends of MHW intensity that are not significant or negligible (Region 1 and Region 2), the values of the two-way information flow both fail to pass the significance test (Figure 2j,l). This also confirms the effectiveness of the information flow analysis in revealing the causality. Therefore, the significant increase trends in MHWs are clearly caused by the increase in atmospheric CO2 concentrations (Figure 1a blue line). However, MHW intensity is above 1 °C in climatology and displays insignificant or even weak negative trends in three regions (Figure 2j-l), contrasting with the slightly increased global mean MHW intensity (Figure 2i). By definition, MHWs are extreme oceanic warming events with large intensity of SST anomalies. MHW intensity is calculated as the mean SST anomalies relative to the climatological daily mean of all the MHWs during specific time intervals, indicating that MHW intensity may also be high during early periods. As MHWs become more frequent and prolonged, the extended MHW days may decrease the averaged intensity of temperature anomalies for MHWs, causing negative trends in the MHW intensity in Region 2 despite the SST increase. Furthermore, Figure 3 shows the spatial pattern of the estimated total change during the 39-year period (1982-2020), based on the linear trend coefficients in annual mean SSTs, and annual MHW days and frequencies in the BOB. The SST displays basin-wide Furthermore, Figure 3 shows the spatial pattern of the estimated total change during the 39-year period (1982-2020), based on the linear trend coefficients in annual mean SSTs, and annual MHW days and frequencies in the BOB. The SST displays basin-wide warming in the BOB (Figure 3a), with a linear trend above 0.4 • C in most regions. Specifically, the warming is slightly high in Region 3 (above 0.5 • C), while it is low in Region 1 (below 0.4 • C). In comparison, there are relatively larger spatial variations in the trends of MHW days and frequencies than the annual mean SST (Figure 3b warming in the BOB (Figure 3a), with a linear trend above 0.4 °C in most regions. Specifically, the warming is slightly high in Region 3 (above 0.5 °C), while it is low in Region 1 (below 0.4 °C ). In comparison, there are relatively larger spatial variations in the trends of MHW days and frequencies than the annual mean SST (Figure 3b,c). Annual MHW days significantly increase east of 85° E, with a trend of above 1.03 days/year, i.e., a total change of 40 days during 1982-2020 that is about 2 times of the 1982-2000 mean value. The MHW frequency significantly increases in Region 3 and south of 12° N in Region 2, but displays an insignificant trend in Region 1. These results indicate that warm-water coral reefs suffer from larger increases in MHW days and frequency in the east and north of the BOB than the west, during 1982-2020. Warm-water reef-building corals need a steady water temperature of 18 °C to 30 °C to live and grow; a surge of 1-2 °C warming above the upper limit of this livable range and sustained over weeks can cause devastating effects on corals. Therefore, coral reefs may be more vulnerable to MHWs occurring at seasons with a background high SST; we further calculated the area-weighted mean MHW frequency in the three regions during boreal winter (December-January-February, DJF), spring (March-April-May, MAM), summer (June-July-August, JJA) and autumn (September-October-November-SON) during 1982-2020 ( Figure 4).
MHW frequency increases significantly in all three regions in DJF, when the seasonal SST (Figure 1e  Warm-water reef-building corals need a steady water temperature of 18 • C to 30 • C to live and grow; a surge of 1-2 • C warming above the upper limit of this livable range and sustained over weeks can cause devastating effects on corals. Therefore, coral reefs may be more vulnerable to MHWs occurring at seasons with a background high SST; we further calculated the area-weighted mean MHW frequency in the three regions during boreal winter (December-January-February, DJF), spring (March-April-May, MAM), summer (June-July-August, JJA) and autumn (September-October-November-SON) during 1982-2020 ( Figure 4).

Interannual Variability of the BOB MHWs
To investigate the interannual variability in warm season MHWs, a 7-year highpass filter is applied to the detrended time series of the March-April-May-June (MAMJ) mean SST and MHW days and frequency in the three regions ( Figure 5). Here, the warmest boreal spring season is selected to highlight the effectiveness of extreme warm water on coral reef bleaching. Indeed, the variations in the MHW days and frequency are highly consistent with the local SST variability at an interannual time scale, with significant correlation coefficients between the local SST and MHW days (frequency) above 0.73 (0.84) in Region 1 and Region 2, and 0.64 (0.77) in Region 3. Specifically, the year ex-

Interannual Variability of the BOB MHWs
To investigate the interannual variability in warm season MHWs, a 7-year high-pass filter is applied to the detrended time series of the March-April-May-June (MAMJ) mean SST and MHW days and frequency in the three regions ( Figure 5). Here, the warmest boreal spring season is selected to highlight the effectiveness of extreme warm water on coral reef bleaching. Indeed, the variations in the MHW days and frequency are highly consistent with the local SST variability at an interannual time scale, with significant correlation coefficients between the local SST and MHW days (frequency) above 0.73 (0.84) in Region 1 and Region 2, and 0.64 (0.77) in Region 3. Specifically, the year experiencing an anomalous high-mean SST also tends to suffer long MHW days and frequent occurrences of MHWs. Furthermore, the years with remarkable anomalous MHWs, characterized by MHW days and frequency exceeding one interannual standard deviation (std), coincide with El Niño years displaying Niño3. 4  El Niño events may associate with the diversity of local atmospheric and oceanic condition changes during different El Niño events. The negative phase of NAO is also suggested to coincide with MHW occurrences in the Indian Ocean [24]. The detailed relationship between the BOB MHWs and typical climate modes has not been systematically studied. Table 1 lists the correlation coefficients between interannual anomalies of MAMJ MHW frequency and typical climate mode indices. Given the dominant role of ENSO in the interannual variability of global climate, the relationship between MHWs and other climate modes is estimated as a partial correlation excluding the effect of ENSO. Here, we use the mean of the Niño3.4 index in December of the developing year and January and February in the decaying year of ENSO, denoted as D(0)JF(1), to highlight the effect of ENSO on the BOB during its decaying phase. For the NAO, previous studies suggest that the NAO in boreal spring can cause a remote effect on the Indian Ocean and East Asia climate in boreal spring and following summer by atmospheric teleconnections [48][49][50][51][52]. Given that NAO peaks in boreal winter, we use the JFMAM(1) mean NAO index to investigate the effect of NAOs. The correlation coefficient between JFMAM(1) NAO index and D(0)JF(1) Niño3.4 index is −0.34, suggesting the necessity of partial correlation analysis (Table 1).   Table 1 lists the correlation coefficients between interannual anomalies of MAMJ MHW frequency and typical climate mode indices. Given the dominant role of ENSO in the interannual variability of global climate, the relationship between MHWs and other climate modes is estimated as a partial correlation excluding the effect of ENSO. Here, we use the mean of the Niño3.4 index in December of the developing year and January and February in the decaying year of ENSO, denoted as D(0)JF(1), to highlight the effect of ENSO on the BOB during its decaying phase. For the NAO, previous studies suggest that the NAO in boreal spring can cause a remote effect on the Indian Ocean and East Asia climate in boreal spring and following summer by atmospheric teleconnections [48][49][50][51][52]. Given that NAO peaks in boreal winter, we use the JFMAM(1) mean NAO index to investigate the effect of NAOs. The correlation coefficient between JFMAM(1) NAO index and D(0)JF(1) Niño3.4 index is -0.34, suggesting the necessity of partial correlation analysis ( Table 1).
As expected, the D(0)JF(1) Niño3.4 index is closely correlated with the MAMJ MHW frequency in all three regions, with significant correlation coefficients of 0.81, 0.65 and 0.53 for Region 1, Region 2 and Region 3, respectively ( Table 1)  MHW frequency with MHW days. The correlation analyses suggest that El Niño and the negative phase of spring NAOs are both important in regulating the interannual variability of MHWs in the BOB, with their relative contribution differing across regions.
The North Indian Ocean SST warming is large at the decaying stage of ENSO, which is also deemed the IOB mode. To explore the seasonal evolution of the relationship between the MHW frequency and Niño3.4 and NAO index, we further calculate the correlation of the D(0)JF(1) Niño3.4 and JFMAM(1) NAO indices with the zonal mean MHW frequency for each calendar month in the three regions ( Figure 6). Specifically, we correlate the D(0)JF(1) Niño3.4 index and JFMAM(1) NAO index with the zonal mean MHW frequency at each calendar month denoted as (1) for 1983-2020, and denoted as (0) for 1982-2019. Note that the correlation between the NAO and MHW frequency is a partial correlation excluding the influence of ENSO. The D(0)JF(1) mean Niño3.4 index is significantly and positively correlated with MHWs from August(0) to September(1) (Figure 6a,b). However, the positive correlation mainly appears from April(1) to autumn(1), coinciding with the second peak warming period of the North Indian Ocean [29][30][31]. These results illustrate the differences of ENSO's effect on MHWs at different latitudes in the BOB, which are important for the prediction of MHWs. When excluding the effect of ENSO, the boreal winter and spring NAOs still display a high correlation with the MHW frequency at certain months in all three regions (Figure 6c,d), with a significant negative correlation mainly appearing from April-May(1) in Region 1 and from January-April(1) in Regions 2 and 3. It is interesting that the influence of the NAO lasts in Region 3 during SON(1) (Figure 6c) when ENSO's effect is relatively weak (Figure 6a). This suggests that the boreal winter and spring NAO may also act as important seasonal predictors for the BOB MHWs.    (Figure 7c, blue line). The correlation coefficient between the reproduced MHW frequency and original MHW frequency slightly increase to 0.85 (Figure 7c). The interannual variability of the two index-regressed time series (std = 0.56, blue line) is also slightly increased by around 4% from that of the Niño3.4 alone-regressed time series (std = 0.54, magenta dotted line), despite that they both underestimate that of the original MHW frequency (std = 0.66). This suggests that the negative phase of the NAO has a limited role in Region 1 MHWs as ENSO is primarily dominant. In Region 2, the influence of ENSO and NAO on the MAMJ MHW frequency weakens as their correlation evidently decreases, compared to that in Region 1 (Figure 7d,f). However, the relative contribution of the NAO compared to the ENSO increases in Region 2, as the percentage increase in the correlation coefficient and the underestimated interannual std of the regressed time series are larger in Region 2 ( Figure 7f) than those in Region 1 (Figure 7c). The effect of the ENSO on MHWs is the weakest in Region 3, compared to the other two regions (Figure 7g), along with minimal interannual std in the Niño3.4 alone-regressed time series (std = 0.31). In contrast, the correlation between the NAO and MHW frequency is the highest in Region 3. There are nearly 30% increases in the correlation coefficients and interannual std of Niño3.4 alone-regression to the joint regression based on Niño3.4 and the NAO. This suggests that the boreal spring NAO is also effective in generating boreal spring MHWs in the BOB, especially in Region 3, possibly through the teleconnection of atmospheric wave trains from the upper level to the lower level [50].

Discussion and Summary
In a previous study [24], the north BOB is noted as the region displaying a high rate of increases in MHW days and frequency during Indian summer monsoon seasons. However, the trend and interannual variability of the BOB MHWs during the climatological warmest season lack investigations, which is the motivation of the present study. Moreover, to highlight the devastating effect of MHWs on coral reef bleaching, we focus the investigation on three regions in which coral reefs are concentrated, with Region 1

Discussion and Summary
In a previous study [24], the north BOB is noted as the region displaying a high rate of increases in MHW days and frequency during Indian summer monsoon seasons. However, the trend and interannual variability of the BOB MHWs during the climatological warmest season lack investigations, which is the motivation of the present study. Moreover, to highlight the devastating effect of MHWs on coral reef bleaching, we focus the investigation on three regions in which coral reefs are concentrated, with Region 1 being located around Sri Lanka, Region 2 in the Andaman Sea and Region 3 in the northeastern BOB. In Figure 5, we note the distinct responses of MHWs during different El Niño years across the regions. Further investigations are needed to explore the underlying processes leading to such differences, with special attention focused on the influences of different ENSO types. The close correlation between the local mean SST and MHW days and frequency suggest that the background monthly mean atmospheric and oceanic conditions caused by climate change are important for the formation of MHWs. The heat budget analyses conducted by, for example, Li and Wang [22] and Saranya et al. [24] are required to fully unravel the physical processes of the BOB MHWs occurring at different locations and during different seasons. The linear trend caused by global warming has not fully emerged for the warm season MHWs in the BOB in recent decades, but once the background temperatures significantly changes, it would fuel more powerful and devastating MHWs. Therefore, it is important to investigate further the future projections of the BOB's MHWs' changes between climate models [8] as they are vital for the adaption of ecosystems and fisheries in the near future.
To summarize, the trend and interannual variability of MHWs in the BOB during 1982-2020 are investigated in this study. In all three regions with concentrated coral reefs, the annual MHW days and frequency significantly increase during 1982-2020 at rates that are larger than that of the global mean. Causality analyses show that the increase trends in annual MHW days and frequency are clearly driven by the rise in atmospheric CO 2 concentrations. However, the linear trend is not significant during the warmest season south of 15 • N due to large interannual variability. In contrast, the coldest boreal winter season witnesses significant and steady increases in MHW days and frequency. These suggest that the most harmful MHWs for the BOB coral reefs are still dominated by natural variability.
The interannual variability of MHWs in all three regions is primarily dominated by ENSO. The boreal winter and spring NAO also modulate their occurrences, especially over the northeastern area in the BOB. The two climate modes synergistically explain about 50~70% of the interannual variances in MHWs in the three regions. Correlation analyses reveal the relationship between MHWs and ENSO and NAO in detail. Specifically, the effectiveness of El Niño in generating MHWs is evident from the boreal autumn of its developing phase to the boreal summer of its decaying phase in Region 1 and Region 2, along with the limited influence of NAOs. However, the effect of El Niño merely appears from April to August of its decaying stage, i.e., the second peak warming peak period in the north of the Indian Ocean. Additionally, the negative phase of the boreal winter-to-spring NAO also exerts a powerful control over the MAMJ MHW frequency. The results suggest that varied climate modes may jointly or separately influence the BOB MHWs during certain seasons and locations, which is important for the prediction of the BOB MHWs. Indeed, when we combine the D(0)JF(1) Niño3.4 index and JFMAM(1) NAO index to build a regression model, it is more effective in reproducing the MHW frequency compared to that of the D(0)JF(1) Niño3.4 index alone, especially in the northeastern area of the BOB.

Conflicts of Interest:
The authors declare no conflict of interest.