A high-resolution (3–8 km) regional oceanic general circulation model is utilized to understand the sea surface temperature (SST) variability of Ningaloo Niño in the southeast Indian Ocean (SEIO). The model reproduces eight Ningaloo Niño events with good fidelity and reveals complicated spatial structures. Mesoscale noises are seen in the warming signature and confirmed by satellite microwave SST data. Model experiments are carried out to quantitatively evaluate the effects of key processes. The results reveal that the surface turbulent heat flux (primarily latent heat flux) is the most important process (contribution > 68%) in driving and damping the SST warming for most events, while the roles of the Indonesian Throughflow (~15%) and local wind forcing are secondary. A suitable air temperature warming is essential to reproducing the reduced surface latent heat loss during the growth of SST warming (~66%), whereas the effect of the increased air humidity is negligibly small (1%). The established SST warming in the mature phase causes increased latent heat loss that initiates the decay of warming. A 20-member ensemble simulation is performed for the 2010/11 super Ningaloo Niño, which confirms the strong influence of ocean internal processes in the redistribution of SST warming signatures. Oceanic eddies can dramatically modulate the magnitudes of local SST warming, particularly in offshore areas where the “signal-to-noise” ratio is low, raising a caution for evaluating the predictability of Ningaloo Niño and its environmental consequences.
A sea surface temperature (SST) warming of ~5 K was observed near the west coast of Australia in the southeast Indian Ocean (SEIO) during the austral summer of 2010–11 (Feng et al. 2013; Pearce and Feng 2013; Kataoka et al. 2014; Marshall et al. 2015). Such strong marine heat waves occurring in the SEIO (Wernberg et al. 2012; Zinke et al. 2014) are named the Ningaloo Niño (Feng et al. 2013) in analog to the Benguela Niño in the Atlantic (Shannon et al. 1986; Richter et al. 2010) and El Niño in the Pacific (Philander 1983). The pronounced warm SST anomaly (SSTA) of Ningaloo Niño causes destructive consequences including severe coral bleaching, massive fish die-off, and dramatic changes of the rainfall pattern in Australia (e.g., Stenseth et al. 2002; Depczynski et al. 2013; Tozuka et al. 2014; Zinke et al. 2014; Narayanasetti et al. 2016; Caputi et al. 2019) and also covaries with other important climate modes such as the Indian Ocean dipole (IOD) and El Niño–Southern Oscillation (ENSO) (Feng et al. 2013; Kataoka et al. 2014; Zhang et al. 2018a; Zhang and Han 2018). Due to its pivotal importance, the Ningaloo Niño is attracting increasing attention of the climate and oceanographic communities. Intensive research efforts have been devoted to understanding its dynamics and predictability (e.g., Feng et al. 2013; Kataoka et al. 2014, 2017; Doi et al. 2013, 2015; Kido et al. 2016; Zhang et al. 2018a).
The existing mechanisms of Ningaloo Niño can be broadly categorized into two groups. The first group highlights the critical role of remote forcing from the Pacific through both the atmospheric bridge and the ocean channel processes (Feng et al. 2013, 2015; Kido et al. 2016; Zhang and Han 2018). As a response to the western Pacific warming under La Niña condition, anomalous low-tropospheric cyclonic winds prevail over the SEIO (Feng et al. 2013; Tozuka et al. 2014). Therefore, the coastal northerly wind anomaly against the climatological mean southerly wind can lead to SST warming by driving coastal downwelling and enhancing the warm-water transport of the poleward Leeuwin Current (LC; Kataoka et al. 2014; Tozuka et al. 2014; Marshall et al. 2015). Meanwhile, downwelling oceanic waves forced by the Pacific easterly wind anomalies can penetrate through the Indonesian Archipelagos and further propagate poleward along the west coast of Australia (Clarke 1991; Clarke and Liu 1994; Meyers 1996; Wijffels and Meyers 2004; Feng et al. 2008), which triggers the SST warming in the SEIO through surface heat convergence and enhancing heat transport of the Indonesian Throughflow (ITF) (e.g., Feng et al. 2013, 2015; Kido et al. 2016).
The second group emphasizes the importance of local air–sea interactions in the SEIO (Marshall et al. 2015; Kataoka et al. 2017, 2018; Zhang et al. 2018a). In the “coastal Bjerknes feedback” (Kataoka et al. 2014; Kido et al. 2016), negative sea level pressure (SLP) anomalies associated with the anomalous cyclonic winds are generated to the west of warm SSTAs, which amplifies the existing SST warming through enhancing coastal downwelling and the LC. The reduced wind speed by the northerly wind anomaly in the presence of mean southerly winds reduces the upward latent heat flux QL and strengthens the warm SSTAs through the wind–evaporation–SST (WES) feedback (Marshall et al. 2015; Zhang et al. 2018a). Kataoka et al. (2017) showed that QL can also damp warm SSTAs if northerly alongshore winds are not well established and proposed a positive feedback effect by the time-varying mixed layer depth. The cloud–radiation–SST feedback, involving the reduced high cloud amount, increased surface solar radiation, and warm SSTA, was also suggested to be essential in regulating SST variability of Ningaloo Niño (Zhang et al. 2018a; Tozuka and Oettli 2018; Li et al. 2019).
Recently studies attempted to evaluate the contributions of different processes in a quantitative manner using oceanic general circulation model experiments (Zhang et al. 2018a; Li et al. 2019). Their results showed that the reduced surface turbulent heat flux QT (QL plus sensible heat flux QS) is likely the leading driver of warm SSTAs of Ningaloo Niño, and the effects of the ITF and wind-driven coastal dynamics are relatively weaker. Zhang et al. (2018a) speculated that in addition to the reduced wind speed, the warming and moistening of the lower troposphere caused by meridional wind advection also contribute to the reduced QT. Unfortunately, in these studies, there were no particular experiments performed to quantify the contributions of several key factors such as radiation, air temperature Ta, and air humidity qa. As a result, there still lacks a systematic, quantitative understanding of the SST variability of Ningaloo Niño. Particularly, the relative importance of Ta, qa, and wind speed in causing QT change has not been evaluated yet.
In addition to large-scale air–sea interaction, the potential impact of ocean internal variability on Ningaloo Niño has never been addressed. Existing studies have revealed pronounced mesoscale eddy variability in the SEIO (Jochum and Murtugudde 2005; Feng et al. 2005, 2007; Jia et al. 2011; Zheng et al. 2018; Mao et al. 2019) and considerable rectification effect on low-frequency oceanic variability (e.g., Trenary and Han 2013; Li and Han 2015; Sérazin et al. 2016). Mesoscale eddies may redistribute SST signatures and reinforce local warming magnitudes, aggregating disastrous impacts on ecosystems if collocating with coral reefs. It is unknown whether ocean internal processes have contributed to the formation of the nearshore warming hot spot of ~5 K during the 2010–11 super Ningaloo Niño event (Pearce and Feng 2013). Similar to the synoptic atmospheric disturbances, ocean internal processes such as mesoscale eddies are largely chaotic in nature and unpredictable for the ocean and climate models. Most models used in the simulation and prediction of Ningaloo Niño did not sufficiently resolve oceanic eddies and precluded their effects. The influence of oceanic eddies on the predictability of Ningaloo Niño and its environmental impact is worthy of investigation.
In this study, we aim to address the unresolved issues proposed above to improve our understanding of Ningaloo Niño. All the key processes controlling SST variability are quantitatively investigated using carefully designed experiments of a high-resolution regional ocean model that can adequately represent mesoscale oceanic eddies. The rest of this paper is structured as follows. In section 2, we describe the datasets and the model utilized for our analysis. In section 3, we evaluate the model’s capability in simulating the SST variability of Ningaloo Niño. In section 4, we devote large efforts to unravel the processes that control SST variability of Ningaloo Niño through analyzing the results of regional ocean model sensitivity experiments. In section 5, we explore the modulation effect of mesoscale oceanic processes on Ningaloo Niño and its predictability, with a particular interest in the 2010/11 super Ningaloo Niño event. In section 6, we present a brief summary and discussion for the primary findings of this research.
2. Data and model
To validate the modeled SST variability, in this study we refer to three SST datasets, including the 2° × 2°, monthly Extended Reconstructed Sea Surface Temperature (ERSST) version 4 data (Huang et al. 2015) from 1993 to 2016, the 1.0° × 1.0° monthly Hadley Centre Sea Ice and Sea Surface Temperature (HadISST) data (Rayner et al. 2003) from 1993 to 2016, and the 0.25° × 0.25° monthly Tropical Rainfall Measuring Mission (TRMM) Microwave Imager (TMI) SST data (Wentz 1997) from 1997 to 2014. TMI is based on high-quality microwave satellite SST measurements that can better resolve the skin temperature anomalies of Ningaloo Niño. It also has higher spatial resolutions than the other two datasets and can better resolve oceanic eddy signatures. Following Kataoka et al. (2014), the Ningaloo Niño index (NNI) is defined as the region-averaged SSTA from 108°E to the coast, 28°–22°S. Here the monthly SSTA is obtained by removing the monthly climatology (seasonal cycle) and the linear trend of the entire time series.
b. Model configurations
The Regional Oceanic Modeling System (ROMS) is adopted to perform simulations and experiments of this study. ROMS is a hydrostatic, primitive equation oceanic model that implements a free surface, horizontal curvilinear coordinate and a terrain-following vertical coordinate system (Haidvogel et al. 2000; Peliz et al. 2003; Budgell 2005). In our recent works, this model was successfully utilized to explore the response of coastal regions to the atmospheric forcing (e.g., Guo et al. 2019). The model is configured to a sector-shaped domain covering most of the SEIO region bordered by the Western Australian coast to the east and Java to the north (Fig. 1a). We adopt spatially varying horizontal resolutions, with the grid size enlarging from ~3 km near the west coast of Australia to ~8 km in the basin interior. As such, the complicated processes in coastal areas including oceanic eddies are better resolved, and the SST variability of Ningaloo Niño is more realistically represented in the model. There are in total of 260 × 280 grid points over the model domain. The minimum water depth is set to be 8 m, and the maximum depth of the model domain reaches ~4500 m (Fig. 1b). The model has 30 terrain-following vertical layers, and the surface and bottom vertical stretching parameters are set to be θs = 6.5 and θb = 0, respectively. The vertical resolution in the top 200 m is refined to better resolve the surface mixed layer processes.
There are three open boundaries over the model domain. The northeastern boundary is fixed between Australian and Java, through which the ITF enters from three major outflow passages of the Indonesian Seas (Lombok, Ombai, and Timor). The southern boundary extends from the southwest cap of Australia to 37.4°S, 98.7°E, so that the poleward LC flows out of the domain. The western boundary is situated at ~96.8°E, where the model domain communicates with the Indian Ocean (IO) interior through the South Java Current (Sprintall et al. 1999), the South Equatorial Current (SEC), and the Subtropical Countercurrent (STCC), from north to south. On these open boundaries, the 1/12°, daily Hybrid Coordinate Ocean Model (HYCOM) plus Navy Coupled Ocean Assimilation global reanalysis data (Bleck and Boudra 1981; Large et al. 1994, 1997) are used as the lateral boundary conditions, and a strong relaxation toward HYCOM reanalysis fields is applied to suppress artificial boundary wave signals. The 0.75° × 0.75° European Centre for Medium-Range Weather Forecasts (ECMWF) ERA-Interim product reanalysis data (Dee et al. 2011) are used as surface forcing fields, involving precipitation, air pressure, surface shortwave and longwave radiations, 2-m air temperature and relative humidity, and 10-m winds. The terms QL and QS, wind stress, and evaporation rate are computed automatically by the model using implementing bulk formula (Fairall et al. 1996; Liu et al. 1979). Taking QL as an example, the bulk formula can be written as
where q indicates the water vapor mixing ratio and W expresses the average value of the wind speed. The interfacial value of the water vapor mixing ratio qs is computed by the saturation mixing ratio for pure water at the sea surface interface temperature Ts:
where qsat is the saturation water vapor mixing ratio at the surface air temperature of Ts, and 0.98 is the empirically derived coefficient. In addition,
where Ta is the air temperature and HR denotes the relative humidity. The other constants used in Eq. (1) such as ρa = 1.29 kg m−3, Le = 2.5 × 106 J kg−1, and Ce = 1.35 × 10−3, denote the density of atmosphere, the latent heat of evaporation, and the transfer coefficient of latent heat, respectively.
Figure 2 compares the ocean current data of HYCOM reanalysis in the ITF exit passages with mooring observations of the International Nusantara Stratification and Transport (INSTANT) program (Sprintall et al. 2004, 2009). A visual comparison suggests that HYCOM reanalysis data can realistically represent the mean structure and changes of the ITF at all the four stations in Lombok (8.4°S, 115.9°E), Ombai (8.5°S, 125.0°E), and Timor (11.4°S, 123.0°E and 11.5°S, 123.0°E) passages during 2004–06, except for an underestimation of ~0.2 m s−1 in the along-strait velocity at the Ombai station. Figure 3 compares the average temperature, salinity, and thermocline depth (isotherm depth of 20°C) of the ITF outflow region (100°–120°E, 20°–10°S) of HYCOM reanalysis with the gridded Argo data (Roemmich and Gilson 2009). In general, HYCOM reanalysis data are able to reproduce the interannual oceanic variations in this region. Particularly, the low-frequency variations of the thermocline, which largely represent the Pacific-origin oceanic waves via the ITF (Sprintall et al. 2003), are faithfully represented by HYCOM reanalysis.
c. Model experiments
To reach a quasi-steady state of the model ocean, a spinup run of 20 years under monthly climatologic surface forcing fields and lateral boundary conditions was first conducted. Then the model was integrated forward under 3-hourly ERA-Interim atmospheric forcing fields and daily lateral boundary conditions from HYCOM reanalysis data from January 1993 to December 2016. The simulation period is largely determined by the availability of HYCOM reanalysis data. This experiment, which contains complete remote and local processes, is supposed to simulate the realistic SST variability of Ningaloo Niño and named the Control (Ctr) run.
We also performed a suite of sensitivity experiments (Table 1) for the same period (January 1993–December 2016) to evaluate the effects of different processes and explore the mechanisms of Ningaloo Niño. Each experiment retains the effects of one or several forcing processes, as reflected in its name. The Exp-B, for example, used monthly climatology in all atmospheric forcing fields but daily lateral boundary conditions as in Ctr. This means that Exp-B retains only the forcing effect from the lateral boundaries, as indicated by the letter B (for boundary) in its name. Note that in our case, oceanic waves propagate into our model domain from the northeastern boundary and propagate out via the western and southern boundaries. Therefore, variability in Exp-B mainly originates from the northeast boundary, and this experiment can roughly measure the ITF effect.
In the other six experiments, different atmospheric variability is exerted in addition to the lateral boundary effect. In Exp-BW, 3-hourly 10-m winds were used, but other forcing fields were still fixed to monthly climatology. Therefore, Exp-BW measures the effects of both boundaries (B) and wind forcing (W). Similarly, Exp-BR uses 3-hourly shortwave and longwave radiation (R) and contains the effects of both the B and R forcings. The remnant four experiments were all designed to diagnose the variability and effect of surface turbulent heat flux QT, which are most important and complicated as indicated by previous studies (e.g., Zhang et al. 2018a; Li et al. 2019). Exp-BWTQ used 3-hourly fields in all the factors affecting QT including winds (W), surface air temperature Ta (T), and relative humidity qa (Q) and was able to fully reproduce the effect of QT on SST. Exp-BT, Exp-BQ, and Exp-BTQ were used to isolate the individual effects of Ta, qa, and their combination on QT and SST variability, respectively. Outputs of all the eight experiments for the period of 1993–2016 are stored as monthly fields.
We also performed a 20-member ensemble simulation for the 2010/11 super Ningaloo Niño event that left serious and irreversible impacts on local marine ecosystems (e.g., Feng et al. 2013; Benthuysen et al. 2014). All the ensemble members (i = 1, 2, … 20) are initiated from 1 January 2010 and terminated on 31 December 2011, using identical surface forcing fields and lateral boundary conditions to those of the Ctr run. The initial condition includes the eastward and northward components of currents, water temperature and salinity, and the sea surface level. For member i, the initial condition is the model output of Ctr run on the ith day of January 2010 but still initiated from 1 January. Due to the nonlinearity of ocean internal instabilities, the small difference in the initial condition can lead to dramatically different patterns of ocean eddies among the members by the peak season of Ningaloo Niño (DJF of 2010/11). As such, the differences among ensemble members quantify the effects of ocean internal processes and are regarded as “noise” while their ensemble mean represents the deterministic, atmosphere-forced variability and is considered as the “signal.” We need to state that the signal and noise here are defined from the oceanic point of view. In fact, the atmospheric forcing also contains noises arising from atmospheric internal instabilities.
d. Model validation for the LC and surface heat fluxes
The experimental strategies described above may cause artificial relaxation effects in the model. Having surface Ta as a boundary condition may tend to keep the simulated SST close to the observed SST, given that SST is close to Ta in reality. Therefore, it is necessary to first verify the model performance in simulating the surface heat fluxes and ocean circulation and confirm whether the model has realistically represented the key processes affecting SST variability. Previous studies showed that the enhanced LC is important for the SST warming of Ningaloo Niño through increased warm water advection to the SEIO (e.g., Feng et al. 2005; Benthuysen et al. 2014; Kataoka et al. 2014; Marshall et al. 2015). Figure S1 in the online supplemental material compares the seasonal climatology of the LC from HYCOM reanalysis data and the Ctr run of ROMS. ROMS can well represent the complicated structure and seasonal variations of the LC, displaying evidently weaker poleward surface currents in January and stronger currents in April and July, respectively. Following Benthuysen et al. (2014), we further calculated the mean surface velocity and the volume transport anomalies of the LC (Fig. S2). As referenced to AVISO geostrophic current and HYCOM reanalysis, ROMS has well reproduced the interannual variations of the LC. The other sensitivity experiments with “wind forcing” included, such as the Exp-BW and the Exp-BWTQ, have also well captured the interannual variations of LC showing the correlations of 0.65 and 0.70 with Ctr, respectively, confirming the dominance of wind forcing and lateral boundary conditions in driving the LC variability. The simulated surface heat fluxes are also consistent with those in ERA-Interim (Fig. S3). The total turbulent heat flux and its primary component (i.e., latent heat flux) are examined. The mean values and seasonal/interannual variations are all in accordance with ERA-Interim data (correlations > 0.86). In addition, the Exp-BWTQ run produces heat fluxes highly consistent with those of the Ctr run (correlations > 0.98), confirming that surface winds, Ta, and qa are the primary factors controlling QT and QL. The fact that the heat fluxes from the Ctr and the Exp-BWTQ are close may not only mean that the radiation anomalies play only a minor role in generating SST variations, but also highlight the dominant roles of the oceanic dynamical response to the surface wind stress and/or the lateral boundary conditions and associated SSTAs in driving the QT variability. The simulation of SST variability is compared with observational data in a comprehensive manner in section 3.
3. The simulated Ningaloo Niño events
The time series of NNI from ROMS Ctr is compared with those from HadISST and ERSSTv4 (Fig. 4). The three NNIs are generally consistent throughout the 1993–2016 period, but the two from ROMS and ERSSTv4 show visibly stronger fluctuations than that from HadISST. For instance, the warming magnitudes in the austral summers of 1999/2000 and 2010/11 are close between ROMS and ERSSTv4 and evidently stronger than that of HadISST. The standard deviation of NNI is 0.508 K in ROMS and 0.519 K in ERSSTv4, whereas it is merely 0.391 K in HadISST. The correlation coefficient between Ctr and ERSST is 0.88, higher than the 0.76 correlation between Ctr and HadISST, although both are significant at 95% confidence. Using the threshold of one standard deviation of the Ctr NNI, eight Ningaloo Niño events were selected during the 1993–2016 period, among which five were “canonical” Ningaloo Niño events (1996/97, 1999/2000, 2010/11, 2011/12, and 2012/13) with their warming peaks occurring in December–February (DJF). The remnant three were “non-canonical” events with their warming peaks occurred in other seasons, such as those in March–May 1998 (1998MAM), March–May 1999 (1999MAM), and September–November 2015 (2015SON). Such noncanonical events have never been reported and compared with canonical events before.
To examine the spatial structure of Ningaloo Niño, the composite maps of peak-season SSTA for the aforementioned eight events are given in Fig. 5. In observations (Figs. 5a,b), the SST warming of Ningaloo Niño reaches the maximal strength (~0.9 K in ERSSTv4 and ~0.7 K in HadISST) near the west coast of Australia and extends northwestward and southwestward to the interior of the IO, and this pattern is well reproduced by ROMS (Fig. 5c). The warming magnitudes of Ctr are close to those of ERSSTv4 in the NNI region (solid rectangle in Fig. 5c). Interestingly, one can see that the warming structure from ROMS, even in the composite of eight events, is much noisier than the observed ones. It shows many regional details and locally enhanced warming hot spots, which are absent in the ERSSTv4 and HadISST composite maps. Since the high-resolution (3–8 km) simulation of ROMS can better resolve mesoscale eddies originating from ocean internal instabilities, it is not surprising to see their strong influence on the warming structure. A natural question arises as to whether the strong eddy influence is realistic or artificially overrepresented by the model.
To answer the above question, we employ microwave SST data of TMI that can better resolve the skin temperature warming signatures with relatively higher spatial and temporal resolutions (0.25° × 0.25° and 3 day, respectively). TMI data covered six Ningaloo Niño events during 1997–2014 (i.e., 1998MAM, 1999MAM, 1999/2000, 2010/11, 2011/12, and 2012/13) and their peak-season SSTA maps are displayed in Fig. 6. As one can see, the SSTA maps from TMI data also show noisy warming structures and mesoscale hot spots in all events, consistent with those of ROMS. Therefore, the mesoscale oceanic eddies indeed strongly modulate the distribution of SST warming signatures and may give rise to warming hot spots that have a strong influence on local marine ecosystems. Since these oceanic eddies are to a large extent unpredictable, they may influence the predictability of Ningaloo Niño, which will be discussed in section 5. Note that the warming magnitudes in TMI, for either large-scale or mesoscale structures, are quite close to those in ROMS. This confirms that ROMS can realistically represent the SST warming of Ningaloo Niño, and HadISST likely tends to underestimate the SST warming amplitude. The good capability of ROMS lends us sufficient confidence to use its sensitivity experiments to explore the mechanisms controlling the SST variability of Ningaloo Niño.
4. Processes controlling SST variability of Ningaloo Niño
In this section, we attempt to clarify the relative importance of different processes controlling SST variability of Ningaloo Niño using sensitivity experiments of ROMS (Table 1). We use the following formula to estimate the contribution of one specific effect to the total NNI anomaly of Ctr run:
where the angle brackets denote the inner product operation, and NNIeffect and NNICtr represent the NNI anomaly caused by an effect (represented by one or a combination of experiments) and simulated by Ctr, respectively. The term C takes into account both the amplitude of NNI variability induced by the effect and its correlation with the total change. First, we examine the time series of NNI derived from sensitivity experiments and compare them with that of Ctr (Fig. 7). As seen in Fig. 7a, Exp-B, Exp-BW, and Exp-BR show evidently weaker NNI variability than Ctr [standard deviation (s.d.) = 0.508 K] because they have excluded some of the forcing processes. Exp-B shows rather weak variations and a standard deviation of 0.174 K. By substituting the NNIs from both Exp-B and Ctr into Eq. (4), we find that the remote forcing effect from the Pacific through the ITF can contribute by ~15% to the total SST variability. When the local wind forcing effect is taken into account, the SST variability in Exp-BW is only slightly enhanced (s.d. = 0.207 K). The local wind effect, including wind-driven ocean dynamical processes (advection, upwelling, and turbulent mixing) and QT changes, can be evaluated as CExp_BW− CExp_B (i.e., 28.6% −15% = 14%). This 14% contains the nonlinear interaction between local wind and ITF. These results are roughly consistent with those derived from other models (e.g., Zhang et al. 2018a; Li et al. 2019). When the forcing effect of surface radiations is added, Exp-BR simulates reduced NNI variability (s.d. = 0.169 K) as compared with Exp-B (s.d. = 0.174 K). This implies that the cloud cover and surface radiations likely provide a weak negative feedback effect (CExp_BR − CExp_B = −1%) to the SST warming, differing from Zhang et al.’s (2018a) analysis that proposed positive feedback by high cloud. Alternatively, our result is in line with Tozuka and Oettli (2018), who showed that the cloud–radiation–SST feedback is positive for Ningaloo Niña (SST cooling) but negative for Ningaloo Niño.
Last but foremost, Exp-BWTQ represents the full effect of surface QT by jointly including of variabilities of winds, Ta, and qa. It reproduces most of the SST variability in NNI with a standard deviation of 0.498 K, which has been quite close to Ctr. To the joint effect of winds, Ta, and qa can be quantified by CExp_BWTQ − CExp_B, which turns out to be ~81%. Note that this 81% contains the effects of both QT and wind-driven ocean dynamics. To better understand the changes of QT, we performed three additional sensitivity experiments (i.e., Exp-BT, Exp-BQ, and Exp-BTQ) to quantify the contributions of different factors (Fig. 7b). This allows for another way to quantify the role of wind forcing effect, that is, using the difference of Exp-BWTQ and Exp-BTQ (CExp_BWTQ − CExp_BTQ = 13%), which is slightly smaller than that estimated using CExp_BW − CExp_B (~14%). This small difference reflects the nonlinear interactions between winds and other forcing processes (ITF, Ta, and qa). If we remove the wind forcing effect (~13%) from the 81% of CExp_BWTQ − CExp_B, the yielded ~68% contribution represents the QT effect induced by Ta and qa changes, which is consistent with the estimation through Exp-BTQ minus Exp-B (69%). Therefore, it is fair to say that the contribution of QT to total SST variability is at least 68%, as wind-induced QT change has been excluded. This provides the first quantitative estimation for its contribution and confirms the dominant role of QT in driving the SST warming of the Ningaloo Niño (Zhang et al. 2018a; Li et al. 2019). In addition, both Exp-BTQ and Exp-BT can capture the majority of the SST variability in Exp-BWTQ, showing standard deviations of 0.445 and 0.441 K, respectively. These experiments and comparisons suggest that surface air temperature Ta acts as the leading contributor to QT and plays a fundamental role in causing warm SSTAs of Ningaloo Niño. Quantified by Exp-BT minus Exp-B, the contribution of Ta to SSTA is ~66% and accounts for the majority of the QT effect. The similarity between Exp-BQ and Exp-B indicates the weak impact of air humidity qa, as confirmed by the much weaker variability in Exp-BQ (s.d. = 0.182 K). Quantified by Exp-BQ minus Exp-B, the contribution of qa is merely ~1%.
The above analysis has clarified the relative contributions of different processes in causing SST variability, and we further look at the event-by-event difference. Figure 8a displays the peak-season NNI of the eight events derived from Ctr, Exp-B, Exp-BW, Exp-BR, Exp-BWTQ, and Exp-BT. In most cases, Ta is of a stronger effect than other factors and generally plays the major role in driving the SST warming. For example, for the 2010/11 super Ningaloo Niño event, Ta contributed ~89% to the total SST warming, which is quite different from Benthuysen et al. (2014), who suggested the first-order importance of the poleward advection induced by the wind-driven LC. In the cases of 1999MAM, 1998MAM, 2011/12, 1999/2000, 2012/13, 1996/97, and 2015SON, contributions of Ta are also dominant, explaining ~82%, ~76%, ~68%, ~65%, ~57%, ~53%, and ~41% of the total warming, respectively. Wind forcing effect, as quantified by Exp-BW minus Exp-B, was important for the 1996/97 (~46%), 2015SON (~25%), 1999/2000 (~24%), 2010/11 (~19%), 1998MAM (~18%), and 2011/12 (~14%) events, but it played a minor role in the case of 1999MAM (~4%) and cast a negative effect on the SST warming in 2012/13 event (−12%). The ITF effect was also considerable in some events, such as 1999MAM, 1999/2000, and 2012/13, in which its contribution was ~31%, ~23%, and ~21%, respectively. It is interesting to note that the effects of Exp-B and Exp-BR generated weak negative SSTAs in the 2010/11 and 2015SON events, implying damping effects in these two events. The three noncanonical events do not show systematic distinctions from the five canonical events in mechanisms. The sample number may be too small to distinguish these two types of events. In addition to the NNI region, a similar analysis is performed to the west box covering the offshore region (see Fig. 5c). The results are broadly similar to those of the NNI region (Fig. 8b and Table S1). Also, Ta acts as the leading contributor to QT and SST changes.
To comprehensively understand the evolution of SSTAs, we further define the developing and decaying stages as the three months prior to and subsequent to the peak season, respectively. Figure 9 shows the composite SSTA maps for the three stages derived from Ctr, Exp-BWTQ, Exp-BR, Exp-BW, and Exp-B. During the developing stage, significant SST warming emerges northwest of Australia in Ctr. The warming is dramatically intensified and extends to higher latitudes in the peak stage. Among sensitivity experiments, Exp-BR, Exp-BW, and Exp-B exhibit evident discrepancies in the amplitude and distribution of SSTAs from the Ctr run in all the three stages, confirming the minor roles played by the ITF, wind forcing, and surface radiations. By contrast, Exp-BWTQ generates quite consistent SSTA signatures with those of the Ctr run in all the three stages, suggesting the dominance of QT in the growth and decay of SST warming.
To achieve more in-depth insights into underlying processes, we further show the composites of net surface heat flux Qnet, QT, QL, and QS in Fig. 10. For all the heat fluxes analyzed here, downward (into the ocean) flux is defined positive. In Ctr (the first column), Qnet has warming effects over the entire SEIO during the developing stage. By the peak stage, cooling effects emerge in the subtropical region south of 30°S and the ITF outflow region in the north, although the latter is below 90% significance level. The cooling effect spreads to the NNI region and is greatly strengthened in the decaying stage. These changes in Qnet can largely explain the rise and fall of SST warming. For comparison, the Qnet composites from Exp-BR are also displayed (the second column), which show few similarities to those of the Ctr run and again suggest the limited role of solar radiation through Qnet. However, in Exp-BR, there are negative Qnet anomalies at low latitudes in the peak stage, strong enough to explain those in the Ctr run. This feature implies a cooling effect by surface solar radiation in that region, suppressing the further growth of SST warming and contributing to the balance between driving and damping effects in the peak stage of warming. As revealed by previous studies, the low-latitude sector of the SEIO has high mean SST and is dominated by the high cloud that provides negative feedback for SST warming by reducing surface solar radiation (Doi et al. 2015; Tozuka and Oettli 2018; Li et al. 2019).
Surface turbulent heat flux QT in Ctr (the third column in Fig. 10) shows strong variations and plays the dominant role in determining Qnet in most areas of the SEIO. The only noticed exception is the low-latitude area in the peak stage. The anomalies of QT and Qnet are opposite there, in supportive of the dominance of high cloud and solar radiation interpreted above. The two components of QT, QL, and QS, are also compared (the fourth and fifth columns). We found that QL anomalies are by far stronger than those of QS and largely consistent with QT anomalies, although QS is also generally in phase with QT and makes a minor contribution. During the developing and peak stages of Ningaloo Niño, both the QT and QL can reach as large as ~7 and ~10 W m−2, respectively. Taking the mean mixed layer depth of H = 30 m in this region (Kataoka et al. 2017), an averaged 8 W m−2 anomaly during developing and peak stages in surface heat flux can cause an SST warming of about ~0.54 K in 3 months from developing to peak stage that is enough to explain the total NNI anomaly, as estimated roughly as Q∆t/cpρ0H, where cp and ρ0 are the specific heat capacity and density of seawater, and ∆t = 3 months. In the bulk formula, QL is mainly determined by the relationship of Ta, qa, and wind speed. We showed the simulated QL from Exp-BWTQ, Exp-BQ, Exp-BT, and Exp-BW in Fig. S4. A visual comparison suggests that QL of the Exp-BWTQ is highly consistent with that of the Ctr run, and Exp-BT shows a larger degree of resemblance to Exp-BWTQ than the Exp-BW and Exp-BQ, particularly in the developing and peak stages. We also checked the composite QT maps (figures not shown), and the results greatly resemble those of QL. Through the analysis presented above, we can reach the conclusion that among all the contributors, Ta is the primary one, and it causes and maintains the SST warming mainly through reducing the upward QL. Different from the developing and peak stages, the negative anomalies of QL (and thus QT) in the decaying stages contain competing contributions from Ta, qa, and winds (Fig. S4). Processes in the decaying stage are likely more complicated and require further inspection.
To track the continuous evolution of Ningaloo Niño, we show the monthly composite NNI and QL anomaly in Fig. 11. Subsequent to the rapid growth from −3 to 0 months, the NNI in the Ctr run goes flat from 0 to +1 month (Fig. 11a), indicating a short quasi-equilibrium state in which the driving and damping effects reach a balance. During this period, the NNI in Exp-BW has been decreasing, implying a cooling effect provided by wind forcing or stronger damping effect that cannot be compensated by the weak warming effect. As seen in Fig. 11b, the QL in the Ctr run displays a rapid growth from the −3 month and reaches the maximum of +12 W m−2 at the 0 month, corresponding to the abrupt warming of SST since the −3 month. By the +1 month, the total QL in Ctr has been negative, for which the contribution of Exp-BW is slightly larger than Exp-BT and Exp-BQ, underscoring the importance of negative QL anomalies in initiating the decay of warm SSTAs. Note that the wind forcing effect here also contains the wind-driven ocean dynamical processes. In the following decaying stage, all of Ta, qa, and winds contribute to the negative QL anomalies and the gradual decreasing SST, as shown in Fig. S4. Figure 11b also shows that preceding the peak warming, especially from the −3 to 0 month, the QL of Ctr generally follows that of Exp-BT, confirming the dominance of Ta effect in the reduced QL loss during the SST warming. In the above analysis, we found varying characteristics in different stages of Ningaloo Niño: Ta acts as the leading driver for reduced QL loss and SST warming in the developing stage; the strong SST warming in the peak stage induces negative QL anomalies that initiate the decay of warming; Ta, qa, and winds all contribute to the increased QL loss and SST cooling in the decaying stage.
The strong forcing effect of Ta is of interest. The importance of Ta warming in promoting NNI warming was first proposed by Zhang et al. (2018a). It is somewhat against our instinct to see a forcing effect of Ta on SST, since in the tropics Ta generally follows the change of the underlying SST. However, during the evolution of Ningaloo Niño, the two are generally in phase, and SST shows a slight time lag to Ta (Fig. S5). This indicates that Ta change is not a passive response to SST, and there may exist two-way interaction between them. We examined the advection process by surface winds in causing Ta changes, as suggested by Zhang et al. (2018a). The correlation coefficient between the Ta tendency and surface wind advection term is as high as ~0.89. However, the advection term is one order larger than that of Ta tendency. It is possible that wind advection is the primary driver of the growth of Ta, but its forcing effect is strongly damped by other processes such as vertical convection and QT into the ocean. A detailed analysis of these processes affecting Ta has been beyond the scope of this paper.
Here we provide preliminary analysis for the surface wind advection term (Fig. S6). It can be seen in this figure that the linear advection terms (anomalous Ta carried by mean winds and anomalous winds along mean Ta gradient) are the primary contributors to Ta warming, while nonlinear advection terms are less important. For the NNI region specifically, the advection by anomalous winds (the third and fourth columns), especially by anomalous meridional winds, is rather effective in causing Ta warming in the developing stage. However, the advection of anomalous temperature by anomalous wind contributed little to the Ta change of the NNI region. These results are broadly consistent with Zhang et al. (2018a). The anomalies and climatology of winds and Ta are shown in Figs. S7 and S8 to aid the understanding of wind advection terms. The prevailing southeasterly mean winds in the SEIO are favorable for spreading Ta warming signatures from coastal region to ocean interior (Fig. S7), while anomalous northeasterly winds along the west coast of Australia will bring a net warming effect over the NNI region through suppressing the cold temperature advection induced by the mean southerly winds (Fig. S8).
We further examine how surface winds initiate the decaying of warm SSTAs. We plot the evolutions of SSTAs from the Ctr run and QL anomalies from the Exp-BW along with anomalous surface wind anomalies in Fig. 12. As a response to the SST warming, cyclonic wind anomalies are seen since the −1 month (Fig. 12b) and reach the peak strength in the 0 and +1 months (Figs. 12c,d). The southern flank of this cyclone involves southeasterly winds that are in the same direction as the mean winds (Fig. S7), and therefore the total wind speed is increased, leading to negative QL anomalies in the southern areas of the SEIO and being likely to partially contribute to decaying of the SST warming. Moreover, Ta and qa anomalies in response to the peak-stage SST warming may also contribute to the increased upward QL loss and play a role in the decay of SST warming. Therefore, the decay of warming in Ningaloo Niño is caused by increased QL loss associated with peak-stage SST warming and its atmospheric signatures.
5. Modulation effects by ocean internal processes
The high-resolution simulations of ROMS have revealed a strong influence of oceanic eddies on the distribution of SST warming signatures of Ningaloo Niño (see section 3). The composite eddy kinetic energy (EKE) for the peak stage of Ningaloo Niño derived from AVISO surface geostrophic current data and Ctr run are shown in Fig. S9, from which one can see an EKE maximum of ~0.1 m2 s−2 in the LC region. With eddy-resolving configurations, ROMS can well reproduce the distribution of EKE and its elevation in the LC region. The relationship between EKE and Ningaloo Niño is rather tight (Fig. S10). The results showed that EKE is elevated by almost one order from O(10−4) to O(10−3) m2 s−2 in the mature stage of Ningaloo Niño. During the negative counterpart of Ningaloo Niño (i.e., the Ningaloo Niña event), the EKE evolution is generally low, staying on the order of O(10−4) m2 s−2 throughout the event. These results imply that mesoscale eddies can be greatly activated by Ningaloo Niño, probably owing to the strong shear instabilities associated with the enhanced LC. These intensified eddy fields can influence SST signatures, as revealed by TMI-derived microwave SST images and ROMS simulations. Our high-resolution model simulation provides an opportunity to evaluate the impacts of eddy “noises” on the predictability of Ningaloo Niño (Doi et al. 2013, 2015) by computing the signal-to-noise ratio (SNR).
Figure 13 shows the spatial distributions of DJF SSTA of the 2010/11 Ningaloo Niño event simulated by all 20 members. One can see that the large-scale warming structures are nearly identical, with strong warming off the west coast of Australia (~2.0 K) extending northwestward in all members. A close inspection also reveals evident differences among them. For example, mesoscale warming patches with SSTA > 2.0 K are also seen in the offshore region, but their location and size are quite random among these members. This further confirms the important modulation effect by ocean eddies on the distribution of SST warming. To quantify this effect and compare it with forced variability, we plot the NNI time series from January 2010 to December 2011 simulated by 20 members (gray lines in Fig. 14a) along with their ensemble mean (black line). The NNI represents the average SST variability of the NNI region, and it turns out to be generally insensitive to the eddy effect, as indicated by a relatively narrow spread range of members. The ensemble SSTA reaches as large as ~1.0 K, one order larger than the O(0.1) K eddy-induced variability. It is likely that the eddy signatures have been largely averaged out in the computation of NNI. However, when the same analysis is conducted for the SSTA time series of the single point (20.73°S, 110.72°E; marked in Fig. 13), the conclusion is notably different. It is found that the spread among members reaches as large as 1.0 K, on the same order of the ensemble mean SSTA at this point (Fig. 14b), suggesting a strong influence of eddies on local SSTA and its predictability.
To achieve a more comprehensive understanding of the eddy effect on SST variability, the SNR of DJF SSTA is computed for the SEIO (Fig. 15). The SNR shows evident spatial variation with a structure resembling Ningaloo Niño (Fig. 5). The SNR is as high as O(10) near the Western Australian coast and gradually decreases toward the ocean interior. The SNR is reduced to <3 in the subtropical portion of the SEIO west of 110°E, probably owing to the weak forced SSTA (Figs. 5 and 13) and strong mesoscale eddy variability in these areas (Jia et al. 2011; Zheng et al. 2018). The pattern shown in Fig. 15 indicates high predictability of local SSTA near coastal areas and degraded predictability in offshore areas. To improve our understanding and prediction skill of Ningaloo Niño and its environmental impact, the effect of oceanic eddies should be further investigated and taken into account in predicting systems. Moreover, the SSTAs induced by eddies can also influence the overlying atmosphere. Mesoscale air–sea interaction (e.g., Small et al. 2008; Putrasahan et al. 2013; Ma et al. 2016) may play a role in the evolution of Ningaloo Niño.
6. Summary and discussion
As an interannually occurring marine heat wave in the SEIO, the Ningaloo Niño has profound impacts on the regional climate and environment. Many mechanisms have been put forward to explain the generation and evolution of Ningaloo Niño, including the remote forcing from the Pacific and the local air–sea interactions. In this study, we use a high-resolution (3–8 km), regional ocean model to investigate the key processes controlling SST variability of Ningaloo Niño. A series of experiments were carefully designed to quantitatively evaluate the relative contributions of different processes to SST warming. The primary findings are summarized as follows.
ROMS can well reproduce the observed Ningaloo Niño events during 1993–2016, including five canonical cases with the warming peak in DJF (1996/97, 1999/2000, 2010/11, 2011/12, and 2012/13) and three noncanonical cases with the warming peak in other seasons (1998MAM, 1999MAM, and 2015SON). The ROMS simulation reveals complicated spatial structures of Ningaloo Niño with evidently mesoscale noises in SST warming signatures, which are absent in ERSSTv4 and HadISST data. These structures are confirmed by satellite microwave SST data of TMI and likely associated with mesoscale eddy effects.
Using ROMS experiments, we quantified the contributions of different processes. The surface turbulent heat flux QT is shown to be the most important driver (>68%) for SST variability, whereas the effect of the Pacific-forced ITF variability is secondary (<15%). The contribution of surface radiations is small. The term QL is one order larger in variability amplitude than QS and dominates QT variability.
To gain insights into the complicated QL variability, four sensitivity experiments were performed to isolate the effects of surface wind speed, Ta, and qa. The results show that surface Ta warming is the most essential factor in causing reduced upward QL during the growth of SST warming. The increased qa is less important and only explains ~1% of the total SST variability, while the local wind forcing effect, including wind-driven ocean dynamics and QT change, exerts a warming effect of ~13%. Further analysis suggests that Ta warming is likely caused by the northerly wind anomalies along the west coast of Australia through suppressing the cold temperature advection due to mean southerly winds.
In the peak season of Ningaloo Niño, cyclonic winds are established over the SEIO in response to SST warming. It involves anomalous southerly wind in the southern sector of the SEIO. Superimposed on the climatological southerly wind, the total wind speed is increased and therefore leads to enhanced upward QL. Besides, Ta and qa anomalies responding to the SST warming in the mature phase cause the negative QL anomalies, exerting a cooling effect on SST and initiating the decay of Ningaloo Niño. During the decaying stage, the contributions of winds, Ta, and qa to SST cooling are comparable.
A 20-member ensemble simulation is performed for the period of 1 January 2010–31 December 2011 to explore the effects of oceanic eddies on the 2010/11 super Ningaloo Niño event. The ensemble members produce radically different distributions of SST warming signatures that are used to quantify the “noise” induced by oceanic eddies and evaluate the SNR. It is found that the NNI is generally insensitive to such noise, but the pointwise SSTA can be strongly modulated by eddies. The SNR is on an order of O(10) near the west coast of Australia and degraded to O(1) in the ocean interior.
Besides the canonical Ningaloo Niño events peaking in DJF, in this study, we also described three noncanonical events with their peaks in other seasons (i.e., 1998MAM, 1999MAM, and 2015SON). These warming surges do not occur in the warmest season and therefore exert less threat to local biosystems. As a result, these events were much less investigated before. Our analysis presented in section 3 did not find evident differences between the two types of events in driving mechanisms. The forcing by QT and Ta plays the dominant role in causing the SST warming over the SEIO region for both types. It can be argued that the sample number in our simulation may be too small to clearly distinguish these two types of events. Further research involving more events may help clarify this issue.
All the analysis of the present study is based on forced ROMS experiments, and this approach has both advantages and disadvantages. These experiments are valuable because they provide unambiguous insights into the key processes controlling SST variability, which cannot be achieved through coupled model experiments. However, the two-way interaction between atmosphere and ocean is not represented in the ocean model. Our analysis from the oceanic point of view leaves an impression that the QT anomaly is the forcing for the warm SSTA, but in reality the anomalies of QT and its controlling factors (Ta, qa, and wind speed) were all within a positive air–sea feedback loop instead of a force-and-response relationship. For example, surface Ta warming is mainly caused by anomalous surface wind advection that is in turn driven by SST warming. Coupled model experiments may be helpful for clarifying these air–sea interaction processes. On the other hand, coupled models are subjected to their own disadvantages. Long-term integration of climate models tends to display an evident bias in the climatology, which brings unrealistic air–sea coupled processes to the simulated Ningaloo Niño (Kido et al. 2016). Coupled models are also subjected to sensitivity to the coupling frequency. Coupling on the frequency of ~1 h requires massive computing resources, while the coupling of ~1 day precludes air–sea feedbacks on diurnal and subdiurnal time scales (e.g., Danabasoglu et al. 2006; Bernie et al. 2008; Li et al. 2013). In addition, we are able to use high resolutions for a regional stand-alone ocean model as in this study, which is quite rewarding in revealing the effects of mesoscale oceanic eddies. However, high-resolution modeling in both ocean and atmospheric components of a climate model will be again quite expensive in computational resources, although mesoscale air–sea interaction is definitely a theme of interest.
In comparing previous studies that suggest that advection plays a major role in some of the most extreme Ningaloo Niño events (e.g., Zhang et al. 2018a), this study has focused on the SST variability, whereas other studies focused on the upper-ocean heat budget or mixed layer heat budget (e.g., Kataoka et al. 2017). A solid heat balance analysis using the model outputs is therefore necessary to better quantify the contribution from different processes, especially the advection effects, in future study.
Another issue worthy of discussion is the relationship with other climate variability modes that may influence the initiation and development of Ningaloo Niño. For example, the SEIO tends to show more anomalous SST warming surges during positive IOD years (Zhang et al. 2018a). A positive IOD is in favor of Ningaloo Niño by casting northerly wind anomalies near the coast of northwest Australia. As a dominant mode of intraseasonal variability in the tropics, the Madden–Julian oscillation (MJO; Madden and Julian 1971) may also influence the Ningaloo Niño (Marshall and Hendon 2014). This is more so during the La Niña years when the phase speed of MJO is decreased so that the MJO can sufficiently interact with the upper ocean (e.g., Wei et al. 2018; Wei and Ren 2019). For instance, when the deep convection center of the MJO passes the Maritime Continent, cyclonic winds will emerge in the SEIO as a response (Gill 1980). The northerly alongshore winds can cause SST warming and trigger Ningaloo Niño through local air–sea feedbacks. More research should be devoted to investigating the impacts of the IOD and MJO on Ningaloo Niño. Another important interannual climate variation in the subtropical south Indian Ocean is the subtropical Indian Ocean dipole (SIOD), with a warm SSTA occurring southeast of Madagascar and a cold SSTA located in the SEIO during its positive phase (Behera and Yamagata 2001; Zhang et al. 2019). The strengthening of the Mascarene high is accompanied by the stronger southeasterly wind over the northeastern SSTA pole and the weaker wind speed over the southwestern SSTA pole. The large-scale anticyclonic wind anomalies drive the cold SSTA over the west coast of Australia, possibly initiating Ningaloo Niña. For example, the decaying stage of positive SIOD in 2017 April–May is the precondition of the 2017 Ningaloo Nina, which develops during September–November (Zhang et al. 2018b). Given the complicated mechanisms and profound impacts of the Ningaloo Niño, there are many related issues left for further exploration, including the modulations of Ningaloo Niño in amplitude, characteristics, and dynamics by decadal and long-term changes of the background climate state in the SEIO (e.g., Du and Xie 2008; Feng et al. 2010, 2015; Li et al. 2017, 2018), possible impacts on the evolution of Ningaloo Niño exerted by the strong variability the SEIO salinity (e.g., Du et al. 2015; Llovel and Lee 2015), and the role of mesoscale air–sea interaction, which are all interesting subjects for future researches.
We thank three anonymous reviewers for providing nice comments. This work is supported by National Natural Science Foundation of China (NSFC) Grants 41776001 and 41806001 and the National Key R&D Program of China (2019YFA0606702). F. Wang is supported by the National Program on Global Change and Air–Sea Interaction (Grant GASI-IPOVAI-01-01). Discussions with Weiqing Han, Ming Feng, Lei Zhang, and Shoichiro Kido are very helpful for this work. ROMS simulations are performed on the Tian-He supercomputer of National Super Computing Center. HYCOM reanalysis data are obtained from the website http://apdrc.soest.hawaii.edu/data/data.php. ERA-Interim data are downloaded from ECMWF interface website https://apps.ecmwf.int/datasets/data/interim-full-daily/levtype=sfc/. ERSST and HadISST are obtained from http://iridl.ldeo.columbia.edu/SOURCES/.NOAA/.NCDC/.ERSST/.version4/ and https://www.metoffice.gov.uk/hadobs/hadisst/, respectively. The TMI SST product was downloaded from the TMI website (http://www.remss.com/missions/tmi). INSTANT mooring records for the ITF are obtained from http://www.marine.csiro.au/~cow074/. AVISO surface geostrophic current data are available at http://marine.copernicus.eu/services-portfolio/access-to-products/. The gridded Argo data are available on the website http://www.argo.ucsd.edu/Gridded_fields.html.