*2.2. Methods*

Zonal and meridional currents (*u*, *v*) observed at 400 m of EC1 were processed to estimate near-inertial kinetic energy. Raw data extracted from the current meters were quality controlled by following the standard procedure for the instrument types [46,47] and then converted to hourly data through subsampling. The minimum current speed of 1.1 cm/s, measured using a rotor current meter (RCM), was treated as a stall and removed from successive processing. Consecutive missing data of less than 6 h were linearly interpolated, and those longer than 6 h were considered as bad data that were excluded from the analysis. Because the depths where the current meters were mounted differed from the nominal depths as the mooring was tilted by the drag due to horizontal currents, a linear interpolation (or extrapolation) was also performed vertically to yield the horizontal currents (*u*, *v*) at 400 m. The NIWs having zonal and meridional components of horizontal oscillations (*uN IW*, *vN IW*) were extracted from (*u*, *v*), by applying a phase-preserving fourth-order Butterworth bandpass filter, with cutoff frequencies of 0.85 *f* and 1.15 *f*; *f* was approximately 0.0505 cph, corresponding to a period of 19.8 h. The amplitude and horizontal kinetic energy of NIWs were computed as p *uN IW*<sup>2</sup> + *vN IW*<sup>2</sup> and *KEN IW*\_*obs* = 0.5*ρ*<sup>0</sup> *uN IW* <sup>2</sup> + *vN IW* 2 , where *ρ*<sup>0</sup> is the reference density (=1025.0 kg/m<sup>3</sup> ). Because the density is not constant to the reference density but varies over time, the near-inertial potential energy as well as kinetic energy needs to be considered to represent the total mechanical energy.

Instead, to quantify the effect of NIW potential energy variations at 400 m, we deduced the time-varying Wentzel–Kramers–Brillouin (WKB) scaling factor as follows [48]:

$$\left[\mathcal{N}(\mathfrak{x},\mathfrak{y},\mathfrak{z},\mathfrak{t})/N\_{\mathbb{O}}\right]^{-1/2}$$

where *<sup>N</sup>* <sup>=</sup> [−(*g*/*ρ*0)/(*dρ*/*dz*)]1/2 is the buoyancy frequency, and *<sup>x</sup>*, *<sup>y</sup>*, *<sup>z</sup>*, *<sup>t</sup>*, *<sup>N</sup>*0, and *<sup>g</sup>* are the zonal, meridional, and vertical coordinates, time, reference buoyancy frequency, and gravity acceleration (set to 9.83 m/s<sup>2</sup> ), respectively. Seasonal variations are clear and dominant in *<sup>N</sup>* at 50 m (mean and standard deviation are 11.9 <sup>×</sup> <sup>10</sup>−<sup>3</sup> <sup>±</sup> 6.3 <sup>×</sup> <sup>10</sup>−<sup>3</sup> rad/s) while decadal and longer-term changes are significant in *N* at 400 m (mean and standard deviation are 1.5 <sup>×</sup> <sup>10</sup>−<sup>3</sup> <sup>±</sup> 2.4 <sup>×</sup> <sup>10</sup>−<sup>4</sup> rad/s) (Figure 4a,c). In this study, *N*<sup>0</sup> was set to 4.6 <sup>×</sup> <sup>10</sup>−<sup>3</sup> rad/s, based on the EN4 data from the upper 500 m in the vicinity of EC1 from 2000 to 2020, by averaging the vertical profiles of buoyancy frequencies (e.g., vertical *N* profiles in summer and winter are compared to show the seasonal variation limited to the upper layers, Figure 3b). Three vertical *N* profiles obtained from the profiling floats within 80 km from EC1 were used to estimate the vertical profiles of WKB scaling factor (Figure 3c,d).

The intraseasonal variation in *KENIW\_obs* was quantified by applying wavelet analysis to *KENIW\_obs*, and intraseasonal band-averaged variance, defined as VARNIW\_obs\_int, was extracted from the wavelet results. In our study, the intraseasonal band was set to 3–100 days, considering previously reported results on the mesoscale eddies in the region, for example, the mean lifetime of mesoscale eddies was estimated to be 95 days [49]. The MATLAB version of the wavelet toolbox [50] was used (http://atoc.colorado.edu/research/wavelets/, accessed on 12 November 2021). Events for high VARNIW\_obs\_int (period high) were defined as the criterion exceeding 1 standard deviation (σ = 0.10) from the mean (µ = 0.05) of the VARNIW\_obs\_int normalized to its maximum (4.0 <sup>×</sup> <sup>10</sup><sup>6</sup> J <sup>2</sup>/m<sup>6</sup> ) over the total period (i.e., VARNIW\_obs\_int exceeds σ + µ ∼ 0.15 for period high after the normalization; red dashed line in Figure 5c). Events for low VARNIW\_obs\_int (period low) were selected to match the total number of event days (505 days), the same as that of period high after sorting VARNIW\_obs\_int in order of magnitude, and the rest of the periods were defined as period neutral (Figure 5). A total of nine events of period high and seven events of period low were identified. Additionally, we defined the annual and decadal (2000s and 2010s) means of VARNIW\_obs\_int for the interannual and decadal variations in NIW kinetic energy.

**Figure 4.** Time series of (**a**) monthly and (**c**) annual mean Buoyancy frequencies at 50 m (grey line, right *y*-axis in (**a**) only) and 400 m (black lines and squares, left *y*-axis) derived from EN4, (**b**) monthly and (**d**) annual mean Wentzel-Kramers-Brillouin (WKB) scale factors at 400 m. Dashed lines in (**a**) and (**b**) indicate temporal means for 400 m. Blue thick solid and thin dashed lines in (**c**) and (**d**) represent decadal means and their standard deviations for 2000s (from 2000 to 2010) and 2010s (from 2010 to 2020), respectively. **Figure 4.** Time series of (**a**) monthly and (**c**) annual mean Buoyancy frequencies at 50 m (grey line, right *y*-axis in (**a**) only) and 400 m (black lines and squares, left *y*-axis) derived from EN4, (**b**) monthly and (**d**) annual mean Wentzel-Kramers-Brillouin (WKB) scale factors at 400 m. Dashed lines in (**a**,**b**) indicate temporal means for 400 m. Blue thick solid and thin dashed lines in (**c**,**d**) represent decadal means and their standard deviations for 2000s (from 2000 to 2010) and 2010s (from 2010 to 2020), respectively.

The intraseasonal variation in KENIW\_obs was quantified by applying wavelet analysis to KENIW\_obs, and intraseasonal band-averaged variance, defined as VARNIW\_obs\_int, was extracted from the wavelet results. In our study, the intraseasonal band was set to 3–100 days, considering previously reported results on the mesoscale eddies in the region, for example, the mean lifetime of mesoscale eddies was estimated to be 95 days [49]. The MATLAB version of the wavelet toolbox [50] was used (http://atoc.colorado.edu/re-

total period (i.e., VARNIW\_obs\_int exceeds σ+ μ~0.15 for period high after the normalization; red dashed line in Figure 5c). Events for low VARNIW\_obs\_int (period low) were selected

J 2 /m<sup>6</sup>

) over the

mean (μ = 0.05) of the VARNIW\_obs\_int normalized to its maximum (4.0 × 10

**Figure 5.** Time series of (**a**) zonal (ݑேூௐ, black line) and meridional (ݒேூௐ, orange line) components of horizontal NIW oscillations, (**b**) KENIW\_obs, (c) VARNIW\_obs\_int at 400 m in 2003 as an example. Red and blue shaded boxes and no shade in (**c**) indicate period high (H2), period low (L1), and period neutral, respectively. Red dashed line in (**c**) indicates σ + μ. **Figure 5.** Time series of (**a**) zonal (*uN IW*, black line) and meridional (*vN IW*, orange line) components of horizontal NIW oscillations, (**b**) *KENIW\_obs*, (c) VARNIW\_obs\_int at 400 m in 2003 as an example. Red and blue shaded boxes and no shade in (**c**) indicate period high (H2), period low (L1), and period neutral, respectively. Red dashed line in (**c**) indicates σ + µ.

to match the total number of event days (505 days), the same as that of period high after sorting VARNIW\_obs\_int in order of magnitude, and the rest of the periods were defined as period neutral (Figure 5). A total of nine events of period high and seven events of period low were identified. Additionally, we defined the annual and decadal (2000s and 2010s) means of VARNIW\_obs\_int for the interannual and decadal variations in NIW kinetic energy.

The MLD is defined as the depth at which the density (ߩ (changed by a given threshold criterion (∆ߩ (relative to that at a reference depth [44,51]. The threshold was calculated from the temperature change (∆T), relative to that at the reference depth, as follows: The MLD is defined as the depth at which the density (*ρ*) changed by a given threshold criterion (∆*ρ*) relative to that at a reference depth [44,51]. The threshold was calculated from the temperature change (∆*T*), relative to that at the reference depth, as follows:

$$
\Delta \rho = \rho \left( T\_{ref} + \Delta T \lrcorner S\_{ref} \lrcorner P\_0 \right) - \rho \left( T\_{ref} \lrcorner S\_{ref} \lrcorner P\_0 \right)
$$

from the EN4 data at the nearest grid to EC1, and ܲ is the pressure at the sea surface (set to zero). The threshold for temperature change and reference depth were set to ∆T = 0.2 ℃ and 10 m, according to Lim et al. [44]. The σ and μ values of the MLD averaged over the period were 23 and 17 m, respectively (maximum of 84 m in January 2011 and minimum of 10.4 m in August 2006). The MLD varied seasonally, with the minimum and maximum values observed in summer (shallow MLD) and winter (deep MLD), respectively, and interannually with less (more) deepening in winters of 2008–2009 and 2014– 2015 (2002–2003, 2010–2011, and 2018–2019) compared to other years (Figure 4). The MLDs estimated by Lim et al. [44] using the data collected within 1° distance from the EC1 mooring site were compared to validate the MLD estimated from the EN4 data; the correlation coefficient between the two MLD time series was 0.76. where *Tre f* and *Sre f* are the temperature and salinity at the reference depth derived from the EN4 data at the nearest grid to EC1, and *P*<sup>0</sup> is the pressure at the sea surface (set to zero). The threshold for temperature change and reference depth were set to ∆*T* = 0.2 ◦C and 10 m, according to Lim et al. [44]. The σ and µ values of the MLD averaged over the period were 23 and 17 m, respectively (maximum of 84 m in January 2011 and minimum of 10.4 m in August 2006). The MLD varied seasonally, with the minimum and maximum values observed in summer (shallow MLD) and winter (deep MLD), respectively, and interannually with less (more) deepening in winters of 2008–2009 and 2014–2015 (2002–2003, 2010–2011, and 2018–2019) compared to other years (Figure 4). The MLDs estimated by Lim et al. [44] using the data collected within 1◦ distance from the EC1 mooring site were compared to validate the MLD estimated from the EN4 data; the correlation coefficient between the two MLD time series was 0.76.

To estimate the inertial response of the upper ocean in the mixed layer to surface wind forcing, a damped slab model (with zonal and meridional momentum equations) [52,53] was applied, using the following equations: To estimate the inertial response of the upper ocean in the mixed layer to surface wind forcing, a damped slab model (with zonal and meridional momentum equations) [52,53] was applied, using the following equations:

$$\frac{\partial u\_{\rm ML}}{\partial t} = f v\_{\rm ML} + \frac{\tau\_{\rm x}}{\rho\_0 H\_{\rm ML}} - r u\_{\rm ML},\\\frac{\partial v\_{\rm ML}}{\partial t} = -f u\_{\rm ML} + \frac{\tau\_{\rm y}}{\rho\_0 H\_{\rm ML}} - r v\_{\rm ML}$$

where *HML*, *r* −1 , *τx*, *τ<sup>y</sup>* , and (*uML*, *vML*) are the MLD, damping time scale, wind stress, and zonal and meridional currents in the mixed layer, respectively. The *r* <sup>−</sup><sup>1</sup> was fixed to 4 days, as considered in previous studies [36,37,54,55], and the time-varying MLD estimated from the EN4 data was used for determining *HML*.

p The amplitude and kinetic energy of the modelled NIWs were calculated as *uML*<sup>2</sup> + *vML*<sup>2</sup> and *KEN IW*\_*model* = 0.5*σθ*<sup>0</sup> *uML* <sup>2</sup> + *vML* 2 , respectively. The intraseasonalband variance of *KENIW\_model*, defined as VARNIW\_model\_int, was calculated in the same manner as that at 400 m by applying wavelet analysis. The rate of wind work (Π) was calculated using the following equation (inner product of surface wind stress and modelled mixed layer currents):

$$
\Pi = \pi\_{\mathfrak{X}\_{NIW}} \mathfrak{u}\_{ML\\_NIW} + \pi\_{\mathfrak{Y}\_{NIW}} \mathfrak{u}\_{ML\\_NIW}
$$

where *τxN IW* and *τyN IW* are the near-inertial band-passed zonal and meridional wind stresses along the meridional line (see Figure 1). Note that the near-inertial band-passed currents (*uML*\_*N IW*, *vML*\_*N IW*) in the mixed layer estimated using the damped slab model (*uML*, *vML*) represent the near-inertial currents in the MLD to estimate Π.

The ray path of NIWs in the spatially varying stratification along 131◦ E (Figure 1) was computed as follows [56]:

$$\frac{2}{3}y^{3/2} = -\frac{1}{\left(2\omega\beta\right)^{1/2}}\int \left(N(y,z,t)dz\right)$$

where *y* is the meridional travel distance, *ω* is the NIW frequency, *N* is the buoyancy frequency estimated using the EN4 data, and *β* is the meridional gradient of *<sup>f</sup>* (<sup>=</sup> *<sup>∂</sup> <sup>f</sup>* /*∂<sup>y</sup>* <sup>∼</sup> 1.8 <sup>×</sup> <sup>10</sup>−11/m/s). The background flow fields were not considered in the calculation of the NIW ray path.

To examine the effect of the mesoscale flow fields on VARNIW\_obs\_int, background conditions were quantified from the satellite altimetry-derived surface geostrophic currents → *U* = (*U*, *V*). The Okubo-Weiss parameter (*α* 2 ), which diagnoses the relative importance of the strain rate and relative vorticity, is defined as [57] follows:

$$\mathfrak{a}^2 = \left(\mathcal{S}\_n^2 + \mathcal{S}\_s^2 - \mathcal{L}^2\right) / 4$$

where *Sn*, *S<sup>s</sup>* and *ζ* are the normal strain *∂U*/*∂x* − *∂V*/*∂y*, shear strain *∂V*/*∂x* + *∂U*/*∂y*, and relative vorticity *∂V*/*∂x* − *∂U*/*∂y*, respectively. When the total strain *S* <sup>2</sup> = *S<sup>n</sup>* <sup>2</sup> + *S<sup>s</sup>* 2 is larger than *ζ* 2 , *α* 2 is positive, showing a saddle shape of the background flow fields. The effective Coriolis frequency was calculated as follows [27,58]:

$$f\_{eff} = \sqrt{\left(f + \zeta/2\right)^2 - S^2/4}$$

To investigate whether there is a statistically significant effect of mesoscale flow fields on VARNIW\_obs\_int, composite analysis was performed by averaging (composite mean) values of Π, *ζ*, *S* 2 , and *α* 2 separately for period high, period neutral, and period low, and compared to address whether they show statistically meaningful difference using Welch's *t*-test (with 95% significance level; *p*-value < 0.05) (Figure 6, Tables 1 and 2). Then, we classified the events into four categories based on the values of *ζ* and *α* 2 (Table 3). The dominance of *S* 2 to *ζ* <sup>2</sup> was determined by the sign of *α* 2 (e.g., when *α* <sup>2</sup> > 0, the *S* 2 dominated *ζ* 2 ).

**Figure 6.** One-month-long composite maps of geostrophic current (grey arrow) superimposed on rate of wind work Π (in colours) during the period high (27 August–26 September 2016, H6), period neutral (1–30 June 2011), and period low (1–31 August 2007, L4). Green square marker represents the EC1. ଶ ଶ **Figure 6.** One-month-long composite maps of geostrophic current (grey arrow) superimposed on rate of wind work Π (in colours) during the period high (27 August–26 September 2016, H6), period neutral (1–30 June 2011), and period low (1–31 August 2007, L4). Green square marker represents the EC1.

**Table 1.** Mean and standard deviation (in bracket) of four condition parameters (Π, ζ, ܵ

during period high. Bold and underlined values are significant with 95% confidence (*p*-value < 0.05). For comparison, composite values for period neutral are shown in the bottom line. **Period** મ **(10−3 W/m<sup>2</sup> )** ા **(/s)** ࢻ **(s/( (× 10−12/s 2 ) Table 1.** Mean and standard deviation (in bracket) of four condition parameters (Π, *ζ*, *S* 2 , and *α* 2 ) during period high. Bold and underlined values are significant with 95% confidence (*p*-value < 0.05). For comparison, composite values for period neutral are shown in the bottom line.

ߙ and ,

)


L3 17–28 June 2005 0.27 (0.24) −0.06*f* (0.01*f*) 0.07*f* (0.01*f*) 2.77 (3.47) L4 14 August–16 September 2007 0.39 (1.08) 0.23*f* (0.02*f*) 0.11*f* (0.03*f*) −80.73 (25.93) **Table 2.** Same as Table 1 but during period low.



**Table 3.** Categories and corresponding events of VARNIW\_obs\_int events using the two condition parameters of *ζ* and *S* 2 for mesoscale fields. The positive and negative anomalies are denoted by plus (+) and minus (−) signs referenced to a zero value for *ζ* and *ζ* 2 for *S* <sup>2</sup> during each event, respectively.
