In this study, 22 years (from 1994 to 2015) of velocity data from three widely-used global ocean reanalysis datasets are utilized. They are the daily Global Ocean Reanalysis and Simulation Version 4 (GLORYS2V4) produced by Mercator Ocean (http://marine.copernicus.eu) (Lellouche et al., 2013), 3-daily Oceanic General Circulation Model for the Earth Simulator (OFES) that driven by National Centers for Environmental Prediction winds (http://apdrc.soest.hawaii.edu/las_ofes) (Sasaki et al., 2008), and 3-hourly Hybrid Coordinate Ocean Model (HYCOM) version 2.2.99DH, which is the Global Ocean Forecasting System 3.1, (http://www.hycom.org) (Chassignet et al., 2009). The GLORYS2V4 has a horizontal resolution of 0.25° and 75 vertical layers, with a layer thickness that increases gradually from 1 m near the surface to approximately 200 m in the deep layer. The OFES has a horizontal resolution of 0.1°, 54 vertical layers, and a vertical resolution of 5 m near the surface to approximately 330 m in the deep layer. While the HYCOM has a horizontal resolution of 1/12°, 41 vertical layers, and a vertical resolution of 2 m near the surface to approximately 1,000 m in the deep layer. Detailed descriptions of these models are available in the above respective links and research papers.

Also included in this study are moored acoustic Doppler current profiler (ADCP) datasets from the Northwestern Pacific Ocean Circulation and Climate Experiment (NPOCE) (http://npoce.org.cn/) and the Tropical Ocean Climate Study (TOCS) (http://www.jamstec.go.jp/) programs (see Figure 1 ) to validate the model outputs. One mooring from NPOCE was deployed at the east of Mindanao Island (8°N, 127.3°E) from 1 December 2010 to 7 December 2012. The moored ADCPs have recorded velocities from the near-surface down to the depth deeper than 1,200 m by using an upward- and a downward-looking 75 kHz ADCPs manufactured by Teledyne RD Instruments (TRDI) that were mounted on the mooring line at about 400 m. More details of this mooring can be found in Zhang et al. (2014). Two moorings, which are provided by TOCS, were deployed at the north of New Guinea Island (2°S, 142°E) from 12 July 1995 to 5 September 1998 and at the equatorial pacific (165°E) from 31 January 1997 to 8 July 2001. Both recorded velocities of upper ~300 m by using upward-looking 150 kHz ADCPs manufactured by TRDI. Detailed information on these two moorings can be obtained in Kuroda (2000) and Kutsuwada and McPhaden (2002), respectively. These three mooring stations were selected to validate the reanalyses and to provide observational evidence for the MUC, NGCUC, and EUC, respectively.

The reliability of the model results was first assessed by comparing them to observations and suggested that, overall, the models reproduce reasonably well the ADCP’s velocity profiles and their variability (see Figure 3 ). However, some notable discrepancies are clearly shown in OFES and HYCOM as compared to the observations, especially in simulating the weakening (strengthening) of the eastward-flowing EUC during the strong El Niño-Southern Oscillation (ENSO) event in 1997 (1998) ( Figures 3I, O ). Moreover, HYCOM seems to overestimate the northward-flowing MUC ( Figure 3M ). In comparison, the velocity profiles of GLORYS2V4 are in better agreement with the observed velocities in all regions.

To further test the validity of the model results in representing the undercurrents, we compared their time series with that of the observations ( Figure 4 ). We used monthly mean velocity data and selected the depths of 650, 200, and 150 m for the moored ADCP’s location deployed at the MUC, NGCUC, and EUC regions, respectively (see Figure 3 ). Based on the data used in these time series, we applied Taylor diagram analysis (Taylor, 2001) to examine how well each model matches the observation using three distinct statistics, i.e., centered RMSD, correlation coefficient, and STD. In the Taylor diagram, the closer the model is to the observation, which is represented by the RMSD, the better the model approximates the observation. This method determines the error (RMSD) that comes from the discrepancies in pattern similarity (correlation) or from the discrepancies in variance (STD). The results of the Taylor diagram clearly show that the GLORYS2V4 is closest to the observations at the MUC and EUC regions with the correlation coefficient of more than 0.6 and 0.8, respectively (see Figures 4B, F ). In the NGCUC region, all models well approximate the observation with the correlation of more than 0.7. The GLORYS2V4 is slightly higher correlation than others. But, in terms of variance, HYCOM is closest to the observation. However, in this region, HYCOM seems to underestimate the eastward-flowing NGCC (see Figure 3N ).

While the OFES and HYCOM datasets have been extensively used in a number of earlier studies to investigate ocean circulation in the WP (e.g., Chiang and Qu, 2013; Wang et al., 2016b; Nan et al., 2019; Wang et al., 2019; Zhang et al., 2020; Zhang et al., 2021), the GLORYS2V4 is less utilized in this region. Nevertheless, our results have shown that among the models GLORYS2V4 best approximates the observations in three different regions and thereby reasonably well reproduces the regional circulation around the WP. Therefore, we believe that the results of this study can provide a better view of the physical processes in the WTP. As such, in the present study, we only utilize the reanalysis data from GLORYS2V4, including the temperature, salinity, and horizontal velocities.

To evaluate the reliability of GLORYS2V4-based LPTM, we applied this analysis to track the trajectory of surface drifting buoys, obtained by the Global Drifter Program (see Elipot et al., 2022 and https://www.aoml.noaa.gov/phod/gdp/hourly_data.php for details). We selected 5 drifters deployed in the WTP to validate the LPTM analysis. The detailed information of the drifters used in this simulation, e.g., initial location, initial date, and tracking time, is provided in Supplementary Table 1 . The surface zonal and meridional velocities from GLORYS2V4 were utilized in this LPTM simulation. Figure 5 shows the trajectories of both drifters (red line) and tracer-tracking results (blue line). In general, the result of LPTM can reasonably approximate the trajectory of the drifters, in agreement with the background surface currents (see Figure 1 ). Therefore, this result confirms the reliability of the LPTM in tracing water parcels in the WTP. The slight discrepancies may result from unresolved submesoscale dynamics, which is beyond the scope of our study, possibly due to spatial and temporal discretization of the model. Nevertheless, we believe that the LPTM is still reliable and the above limitation will not significantly change the final destination of the tracer.

Before going into the results of the LPTM, it is beneficial to identify the mean structure of the currents that are possibly linked to the NGCUC. Figure 6 displays the vertical profiles of velocity overlaid with the contours of potential density at several sections. The results presented in this figure are based on model outputs from 1994 to 2015, where the first-, second-, and third-row panels are time average and composite during El Niño and La Niña phases, respectively. The sections were selected to show the typical structure of NGCUC and to preliminary identify its potential destinations. There are 180°E, 140°E, 125°E, 8°N, and 1°S sections to show the equatorial eastward currents (e.g., NECC, EUC, NSCC, SSCC, and EDJ), the NGCUC, the ITF, the MUC, and the HTF, respectively. In general, the vertical distribution of mean velocities from the model results shows consistency with previous studies (e.g., Qu et al., 2012; Wang et al., 2019; Delpech et al., 2020; Li et al., 2020; Zhang et al., 2020), displaying the typical formation of zonal or meridional velocities along each section.

The mean state of zonal velocity along 140°E shows that NGCUC ranges from nearly surface down to almost 800 m (or 27.2σ _θ ) and extends almost 200 km from the coast ( Figure 6B ). Please note that the wider range of westward flow near the surface is due to merging with SEC, instead of the wider NGCUC (Wang et al., 2019). The maximum of the mean westward velocity of approximately 68 ± 16 cm s^-1 is found around the depth of 163 m (~25σ _θ ), which represents the NGCUC core.

In 180°E section ( Figure 6A ), a series of alternating zonal currents in the equatorial region is clearly shown, generally consistent with the results of previous studies (e.g., Delpech et al., 2020; Li et al., 2020). The result shows five eastward (i.e., NECC, EUC, NSCC, SSCC, and EDJ) and two westward (i.e., SEC and EIC) currents. The eastward currents shown in this figure are suspected to be the escape routes of the NGCUC waters into the central Pacific. Note that there might be some currents that do not appear in the mean state. For example, the westward core of EDJ (Delpech et al., 2020) and the westward-flowing NESC ( Figure 6K , Li et al., 2020).

Multiple cores of ITF are shown along 125°E section ( Figure 6C ). The main core, which is centered at approximately 4.75°N, is a surface-intensified westward velocity (above 27σ _θ ) with the maximum mean velocity reaching 60 cm s^-1. This core is mainly sourced from the MC that turns westward at the southern tip of Mindanao Island (Feng et al., 2018). Another core is much weaker (~3 cm s^-1) and is centered at approximately 3.25°N and 27.2σ _θ . The source of this core is unknown yet but could be from the equatorial Pacific via NESC (Li et al., 2021) or from the South Pacific via NGCUC. The lower core of ITF along 125°E section will be further investigated in Section 3.3.

Figure 6E shows the vertical profile of southward-flowing HTF. In the mean state, this flow is quite strong with the maximum mean velocity reaching 68 cm s^-1 near the surface. Along 1°S section, the HTF core is mainly centered at the western part of the channel (~128.5°E) along the east coast of Halmahera Island. The result in Figure 6J suggests that the velocity core of the HTF is revealed deepening during El Niño phase. A detailed analysis of this event will be described in section 3.5.

Section 8°N shows a core of northward-flowing MUC with the maximum mean velocity reaching 4.5 cm s^-1 at the depth of around 560 m ( Figure 6D ). This core extends approximately 120 km offshore and is found below 400 m (or below ~26.75σ _θ ). Previous studies have confirmed that this current is a quasi-permanent current with strong intraseasonal variability, which is closely related to SE-1 and SE-2 activities (Azminuddin et al., 2022).

The destination of NGCUC is deduced from the trajectories of tracers released from the NGCUC region using LPTM. Figures 7A–F display the resulting trajectories that represent the pathways of tracers in all selected density layers. The color shading denotes the mean probability distribution of NGCUC water parcels at each grid cell (0.25°×0.25°). The probability was estimated as the ratio between the number of tracers passing through grid cells and the total number of tracers released. The grid cells with higher probability are the region where the tracers pass more frequently. This interpretation can identify the destination of NGCUC including predicting its water mass distribution. The results show that the dissemination of NGCUC tracers varies markedly with depth, but in general, most tracers are distributed to the east, which is consistent with previous studies (Wang et al., 2016b; Li et al., 2020; Zhang et al., 2020; Li et al., 2021).

As described in Section 2.2, some predefined sections were selected to quantitatively classify the destination of each tracer (see Figure 2 ). The classification of the NGCUC tracers’ fate from all isopycnal layers is summarized in Figures 7G-I . Within two years’ travel time of tracer tracking, approximately 50.48%, 0.62%, and 30.84% of tracers are classified as major, minor, and transient destinations, respectively. Note that the total classified tracers into these destinations explain only 81.94% of the total tracers that were released daily from 1994 to 2013. The classified tracers are those that reach a predefined section, meaning the remaining unclassified portions (18.06%) do not pass through any predefined sections within the given travel time, possibly due to being slowed down by turbulence or the effect of a random-walk component in Eq. 1 that could suddenly change the direction of tracers’ motion before reaching any destinations.

The tracers classified as the transient destinations are those that are trapped in the Pacific Western Boundary Eddies (hereafter called the WBE, defined here as the HE and SE-2 regions) or the Bismarck Sea (BS). The present study considers both WBE and BS as transient destinations (or transit locations) rather than major destinations of NGCUC, in which the tracers need more travel time to escape from there and be redistributed into one of the major or minor destinations. The results show that approximately 22.76% (8.08%) of total tracers are trapped in the WBE (BS).

To confirm this assumption and to further identify the fate of tracers trapped in those transient destinations, we released tracers in both WBE and BS and applied the same scenario of LPTM (see Figures 7G-II-III ). It implies that the portion of NGCUC tracers trapped in the WBE and the BS will be distributed according to these respective ratios. The maps of mean probability distribution of tracers released at these regions can be seen in Supplementary Figures 1 and 2 . The result reveals an overall decrease (increase) in the tracers distributed to the eastern destinations for the WBE (BS)-originated tracers as compared to the NGCUC-originated tracers. Notably, the tracers distributed to the MUC account for more than one-fourth of all tracers released from the WBE ( Figure 7G-III ) – a much bigger slice of the destination pie than the tracers directly released from the NGCUC region (i.e., 1.94%) ( Figures 7G-I ). Consequently, the total ratio of NGCUC’s tracers distributed northward through the MUC increases up to 11.3% becoming 13.24% ( Figure 7G -IV). A large portion of the WBE-originated tracers distributed to the MUC region is also clearly shown in Supplementary Figures 2D–F . This implies a critical role of eddy in trapping and transporting the NGCUC waters, and in connecting the NGCUC to the MUC, which is typically disconnected. A more detailed analysis of the WBE will be provided in Section 4.

By considering all NGCUC tracers directly originated from the NGCUC region, and the tracers temporarily trapped in the WBE and BS, the present study suggests the total destinations’ ratio of the NGCUC tracers as summarized in Figure 7G-IV and Table 1 . The results confirm that a large amount of the NGCUC waters is distributed to the central Pacific through the EUC, NECC, NSCC, SSCC, and EDJ. It accounts for 81.65% of all tracers released. Most of the remaining portions are distributed westward into the MUC, ITF, and HTF. This result may reflect the distribution’s ratio of the total NGCUC transport.

By means of LPTM, various destinations of NGCUC are observed, but only eight major destinations are further analyzed in the present study. There are NECC (12.3%), EUC (35.26%), NSCC (13.33%), SSCC (8.85%), EDJ (11.49%), MUC (13.24%), ITF (3.47%), and HTF (0.86%). These are suggested as the major destinations considering the relatively high percentage (>0.5%) and the persistency of tracers passing through their territory or tracks. In addition, it is also observed that a few portions of tracers move toward the Maluku Sea (MS) (0.43%) and the Solomon Sea (SS) (0.42%), and the KC (0.34%). The tracers distributed into these minor destinations (i.e., SS, MS, and KC) explain approximately 1.2% of the total tracers released. The finding of these minor destinations, especially the KC and the SS, is surprising and beyond what we expected at the beginning of this study. Note that the upper NGCUC waters can be advected back to the SS by the southeastward-flowing NGCC. We suggest that some tracers, especially above 25σ _θ , can be advected northwestward to the KC after joining the NECC to the east and then the NEC to the west possibly due to the strong NEC phase (see Figures 5 , 7A, B ). This implies a possible connection between the NGCUC and the KC. However, further investigations with observational evidence are required to confirm this finding.

The results of LPTM show that different isopycnal layers can have disparate distributions of NGCUC waters. Along the isopycnal layers of 23 and 25σ _θ , only about 75.77 and 82.57% of total tracers, respectively, are successfully classified as major destinations. Most of the remaining portions are instead being trapped in the HE and the BS. A similar situation may also apply for the deeper layer, but being trapped in the SE-3 and the BS. Please note that the percentage of classified tracers is not only affected by the eddy’s and the BS’s trapping effects but also by the velocity of the flow itself. At the deeper layer, the current’s velocity becomes slower and thus the tracers need more time to reach any destinations’ section. For this reason, the total classified tracers are getting lower at the deeper layer (see Table 1 ). The percentage of classified tracers may increase with longer travel time, but we believe that it will not significantly change the ratio of the major destinations.

In the probability map ( Figures 7A–F ), the spreading of tracers is clear. The results show that the majority of NGCUC’s water parcels turn eastward to the central Pacific through several eastward-flowing currents. However, their portions vary by isopycnal layers. Along the uppermost isopycnal surface (i.e., 23σ _θ ), more than 67% (19%) of total classified tracers are distributed toward the east into NECC (EUC) around 5°N (equator). In comparison, at the isopycnal layers of 25 and 26σ _θ , most of the tracers (82.5% and 95.6%, respectively) turn eastward along the equator representing the core of EUC. At the deeper layer (26.75σ _θ ), this prominent eastward core disparts into northern and southern branches representing NSCC (54.6%) and SSCC (30.5%), respectively. The probability map shows that a core of equatorial eastward flow re-emerges at the isopycnal layer of 27σ _θ representing EDJ (50.7%), while the NSCC and SSCC branches are vastly lessened. These three eastern destinations keep appearing at the deepest isopycnal layer (27.2σ _θ ) but with a shorter zonal extent and their cores seem to merge and are unclearly distinguished. Please notice that the vertical and latitudinal ranges of these eastern destinations are consistent with the result of the mean zonal velocity distribution shown in Figure 6A . Their cores are well represented by the tracer distribution, implying the reliability of the LPTM result.

The portions of tracers that are distributed to the western destinations (i.e., MUC, ITF, and HTF) also vary by isopycnal layers. Up to 13.24% of total classified tracers are distributed poleward along the Philippine coast joining MUC. This portion is evident at the isopycnal layers below 26.75σ _θ (see Figures 7D–F ) and the tracers may reach up to 15°N. This probability distribution is consistent with the poleward extent of AAIW (up to 15°N along the Philippine coast), which originates from the South Pacific via NGCUC, reported by Qu and Lindstrom (2004).

Around 3.47% of tracers are distributed to the ITF entrance (Celebes Sea) mainly through the south of Talaud Islands. These tracers are evident in most layers but mainly distributed at the uppermost (23σ _θ ) and lowermost (27.2σ _θ ) isopycnal layers which explain about 11.2% and 5.2% of total classified tracers, respectively. This result is consistent with the distribution of mean westward flow along 125°E section that also shows the westward cores around both 23σ _θ and 27.2σ _θ isopycnal layers (see Figure 6C ). The LPTM results confirm that the southern westward flows (2.5° – 4°N) are mainly sourced from the NGCUC, while that at the northern (4° – 5.5°N) is obviously mainly sourced from the MC.

Other portions of NGCUC tracers also enter the Indonesian Sea, but through the eastern passage via Jailolo Strait (i.e., HTF). The number of tracers crossing HTF’s section was only 0.86% from the total classified tracers. Although the percentage is relatively small, the results from LPTM confirm that the NGCUC tracers persistently pass through this route, especially at the upper layers. After entering the Jailolo Strait, the tracers flowed toward the eastern part of ITF route, including the Halmahera Sea, Seram Sea, Lifamatola Passage, and Banda Sea before finally arriving in the Indian Ocean through Ombai Strait or Timor Passage. South Pacific water advected into these Indonesian Seas through the HTF has been previously confirmed (e.g., Gordon et al., 2003), which is consistent with our results.

In Section 3.3, various destinations of the NGCUC have been identified. Nevertheless, the ratio of the destination may change over time. To identify the major signals that significantly control the spatiotemporal variation of the NGCUC destination, we applied the power spectral density function for the transport of the NGCUC and its major destinations. In most regions, the power spectrum results show peaks at seasonal and interannual periods (see Supplementary Figure 3 ). Further analyses using LPTM also confirm that the ratio of NGCUC destinations evidently changes with both timescales. Therefore, this study will focus mainly on seasonal and interannual variations.

Although LPTM has been able to estimate the spatiotemporal change of NGCUC’s destinations, the real distribution of NGCUC water mass at each destination has remained to be clarified. We then further analyze the water mass characteristics at each region by using θ -S diagram. Figures 12A, B display θ -S diagrams at the NGCUC region and its major destinations, which are estimated from model outputs and averaged from 1994 to 2015. To estimate the θ -S diagram, a point, where tracers pass more frequently, is selected as a proxy for each destination (i.e., 1.5°S, 140°E [NGCUC]; 8°N, 127.5°E [MC/MUC]; 1°S, 128.75°E [HTF]; 3.25°N, 125°E [ITF]; 3.5°S, 148.5°E [BS]; equator, 160°E [EUC/EDJ], see Figure 2 ).

Based on θ -S diagram analysis, the upper NGCUC water is characterized by a salinity maximum core (>35.1 psu) centered at 25σ _θ representing the SPTW that originates from the subtropical South Pacific (Fine et al., 1994; Qu et al., 1999). At the intermediate layer, NGCUC carries AAIW, which is characterized as salinity minimum water (34.55 psu) at 27.2σ _θ , originating from the Antarctic convergence zone (Qu and Lindstrom, 2004). The signatures of SPTW and AAIW, which are carried by NGCUC, can be found in most of the selected θ -S diagram points but with various transformations. Along 25σ _θ , the SPTW carried by the NGCUC is transformed into fresher water at BS and EUC regions. While at the other regions (i.e., HTF, ITF, and MC), the upper-layer water is mainly characterized by the fresher water NPTW, instead of the SPTW. It implies that upper-layer waters at the Pacific equatorial western boundary are mainly sourced from the MC that carries the NPTW.

At the MC/MUC region (8°N, 127.5°E), the θ -S diagram is characterized by the North Pacific waters, which is obvious, including the North Pacific Tropical Surface Water (NPTSW, a surface salinity minimum), the North Pacific Tropical Water (NPTW, a subsurface salinity maximum in the thermocline), and the North Pacific Intermediate Water (NPIW, a salinity minimum in the subthermocline) (Grey line in Figure 12B ). However, below 27σ _θ the South Pacific water, i.e., AAIW, is confirmed, which is consistent with previous studies (e.g., Bingham and Lukas, 1994; Fine et al., 1994; Wang et al., 2016a). The results in Figures 6D and 12B suggest that in this region the fresher water NPIW is carried by the southward-flowing MC and the saltier water AAIW is transported by the northward-flowing MUC below 26.75σ _θ . The possible mechanism of how the South Pacific waters are transported to the Philippine coast has been suggested in our previous study (Azminuddin et al., 2022), which is carried by the NGCUC and the SE-2. This result is also well represented by the LPTM results in Section 3.3, implying the reliability of the particle tracking method.

Figures 12C–H show the θ -S diagrams with climatologically monthly mean (color θ -S diagrams) in all locations. The seasonal change of the water mass properties is clearly found in some regions, especially in the HTF region ( Figure 12E ). The results show that at the upper 24.5σ _θ , the HTF water transforms from relatively warm and salty water in summer to colder and fresher water in winter. This implies that during summer when the NGCUC is strengthened, the upper layer of HTF is transformed into saltier water possibly carrying the SPTW from NGCUC. While in winter, the upper layer of HTF is mainly sourced from the NPTW, which is characterized by relatively cold and fresh water. Note that this seasonal change is altered around the isopycnal layers of 25 – 25.7σ _θ . The results of LPTM further confirm that at the isopycnal layer of 23σ _θ the tracers are distributed more to the HTF in summer than in winter, but the case is reversed at the isopycnal layer of 25σ _θ (see first- and second-row panels in Figure 9 ). This result suggests that seasonal changes in the NGCUC strength can substantially magnify the seasonal water mass transformation of the HTF.

The interannual change of the θ -S diagrams at each region is also estimated as shown in Figures 12I–N . The black, red, and blue lines denote the time average and composite during El Niño, and La Niña, respectively, from 1994 to 2015. The AAIW at the MUC region is clearly defined in the La Niña phase. But the clarity is slightly lessened, which is fresher during El Niño years possibly mixing with the fresher water of NPIW. Figure 6I confirms that during El Niño the MC is strengthened and deepened, shifting the MUC further offshore. This event may transport more North Pacific waters to the MUC region and transform the intermediate water mass into becoming NPIW.

Figure 12K displays a noticeable water mass transformation of the HTF associated with the ENSO cycle, in which the upper HTF water is transformed into colder and fresher (warmer and saltier) water in the El Niño (La Niña) phase. This result implies that in addition to the NPTW and SPTW, the upper HTF can be sourced from another region. During La Niña (stronger NGCUC), more SPTW waters enter the HTF resulting in warmer and saltier water. However, during El Niño the source of fresh and cold waters that lead to a unique water mass remains uncertain. It can be either from NPTW or Indonesian Seas (e.g., Banda Sea). Koch-Larrouy et al. (2006) confirmed that the Banda Sea is characterized by a relatively constant salinity (around 34.58 psu) below 20°C. This is consistent with the upper HTF water during El Niño phase (red line in Figure 12K ). Moreover, the velocity core of the HTF is revealed deepening during El Niño phase (see Figure 6J ), and its velocity upper 25.5σ _θ reverses becoming a northward current, flowing from the Indonesian Sea to the Pacific Ocean. So, we suggest that during El Niño the upper HTF water mainly comes from the Banda Sea.