Spatiotemporal variability of the subduction rate in the North Pacific from 2005 to 2012 is examined based on the Argo observational data. The subduction rate in the subtropical North Pacific varies significantly from year to year between 25 and 50 Sverdrups (Sv; 1 Sv ≡ 106 m3 s−1), and it is well correlated with the Pacific decadal oscillation. The temporal change of the subduction rate is largely determined by that of the late winter mixed layer depth through the lateral induction term. The increase (decrease) in the subduction rate in the subtropical mode water areas accompanies densification (lightening) of the mode density class of the subducted water. The subduction rate variability in the central mode water and eastern subtropical mode water regions is anticorrelated as found in the previous study using the output from an ocean GCM. The subduction rate in the central mode water density range changes dramatically, which is very large in 2005 and 2010 but almost disappears in 2009. The subduction rate variability in the western subtropical mode water regions seems to be correlated with the Pacific decadal oscillation with a lag of a few years.
The subduction process in which water is transferred from the surface mixed layer into the permanent pycnocline is a key process to understand the formation of the permanent pycnocline and subsurface water masses such as mode waters. The upper ocean is primarily ventilated through the subduction process (Williams 2001). For this reason, subduction rates are routinely calculated to quantify the ventilation process, to elucidate the link between surface and subsurface water properties, and so forth.
The mode waters are formed as a deep mixed layer in late winter due to convective mixing and are left in the subsurface after being capped by the seasonal pycnocline in spring (Hanawa and Talley 2001; Oka and Qiu 2012). The mode waters are generally formed on the warm side of a strong front and appear as a subsurface thermostad or pycnostad, displaying vertical homogeneity with respect to their water properties. Part of the formed mode waters enters the permanent pycnocline (i.e., subducts) and spreads over a wide area through advection. In the open North Pacific, there are at least four mode waters: Western Subtropical Mode Water (WSTMW; Masuzawa 1969), Eastern Subtropical Mode Water (ESTMW; Hautala and Roemmich 1998), Central Mode Water (CMW; Nakamura 1996; Suga et al. 1997) and Transition Region Mode Water (TRMW; Saito et al. 2007). CMW is said to have lighter (L-CMW) and denser (D-CMW) variants in association with frontal structure in its formation region (Oka and Suga 2005).
Oka and Qiu (2012) reviewed recent progress in North Pacific mode water research and summarized the latest knowledge on mode waters such as their formation areas, water properties, and advection routes. WSTMW is formed south of Kuroshio and the Kuroshio Extension Front (KEF) between 132°E and near the date line. L-CMW (D-CMW and TRMW) is formed in a latitude band of 33°–39°N (39°–43°N) elongated from 142°E to 160°W. D-CMW (TRMW) is predominantly formed in the eastern (western) part of the formation domain. The formation of these mode waters is related to the Kuroshio Bifurcation Front (KBF) and Subarctic Front (SAF). Oka et al. (2011) analyzed profiling float data obtained during 2003–08 and showed the density classes of WSTMW, L-CMW, D-CMW, and TRMW are 24.2–25.6σθ, 25.4–26.3σθ, 26.1–26.5σθ, and 26.4–26.6σθ, respectively. As schematically summarized in Fig. 2 of Oka and Qiu (2012), these mode waters are advected on isopycnal surfaces. WSTMW is thought to be advected southwestward mainly by recirculation gyre flow. L-CMW is advected southeastward and then southwestward along large-scale subtropical gyre flow. D-CMW and TRMW is advected along outer paths of the gyre periphery. ESTMW is formed as a relatively deep late winter mixed layer, associated with the subtropical/subpolar water mass boundary near 25°–30°N, 135°–140°W, which is advected southwestward after subduction (Hautala and Roemmich 1998).
The subduction rate has been calculated using ocean climatology based on historical observation data. Huang and Qiu (1994) estimated the subduction rate in Lagrangian coordinates using Levitus climatology data (Levitus 1982) and Hellerman and Rosenstein wind stress data (Hellerman and Rosenstein 1983). The total subduction rate in the North Pacific subtropical gyre was estimated to be 33.1 Sverdrups (Sv; 1 Sv ≡ 106 m3 s−1). Suga et al. (2008) calculated the subduction rate using the same method with Huang and Qiu (1994) but using HydroBase (Macdonald et al. 2001) and wind stress climatology from ECMWF ERA-40 (Uppala et al. 2005). They estimated the total subduction rate as 33.6 Sv (37.8 Sv) when the North Pacific subtropical gyre is defined as the region of the negative wind stress curl north of 15°N (as the region of the positive Sverdrup streamfunction north of 15°N). They used various wind stress climatologies to check the sensitivity of their estimate and found that the Ekman pumping velocity from the various climatologies do not differ much in terms of their effects on the annual subduction rate. The total subduction rates estimated by Huang and Qiu (1994) and Suga et al. (2008) using different ocean climatologies and various wind stress climatologies are all quite similar. An ocean general circulation model (OGCM) has been also used to examine the climatological aspect of the subduction rate. Qu et al. (2002) used an output from an eddy-permitting (¼° in both latitude and longitude) OGCM to evaluate the subduction rate in the model. The total subduction rate in the North Pacific (10°–50°N) was estimated as 61.6 Sv, which is much larger than those derived from the climatologies (Huang and Qiu 1994; Qiu and Huang 1995; Suga et al. 2008). They argued that the larger subduction rate may be attributed not only to mixed layer parameterization and model biases (e.g., a warm bias of the Kuroshio and its extension of their model, a bias common to models that restore SST back to the observed climatology) but also to the higher spatial resolution of the model. Smoothing employed to construct the observational climatologies is generally heavy, which weakens the mixed layer depth (MLD) gradient and reduces the lateral induction substantially.
Time-dependent features of the subduction rate have recently been studied based on OGCM outputs and assimilation products. For example, Liu and Huang (2012) used the Simple Ocean Data Assimilation (SODA) outputs to estimate the global subduction rate and its variability over the period from 1959 to 2006 and found substantial variability of the subduction rate on interannual and decadal time scales, which is predominantly caused by MLD variability through lateral induction. The total subduction rate in the North Pacific averaged over the period is estimated as 82.4 Sv. The larger value compared with previous studies using climatologies is partly due to the larger domain of integration. They also suggested that the usage of climatologically mean fields rather than year-to-year varying fields may greatly underestimate the annual subduction rate because of the nonlinearity associated with the lateral induction term. Qu and Chen (2009) examined the subduction rate based on the OGCM for the Earth Simulator (OFES; Sasaki et al. 2004) output with a horizontal resolution of 0.1° in both latitude and longitude and obtained a total subduction rate in the North Pacific (20°–50°N) of 49.3 Sv. They found large decadal variability in the subduction rate, with these variations corresponding well with the Pacific decadal oscillation (PDO; Mantua et al. 1997). Chen et al. (2010) further examined this decadal variability in subduction rates in the North Pacific by dividing it into western and eastern components and found that the variability of the two components is anticorrelated. That is, the subduction rate in the western (eastern) North Pacific is large (small) during the warm phase of PDO, which is primarily caused by deeper (shallower) MLD in the west (east) during that phase of PDO.
Despite advances in characterizing temporal variability in subduction rates using OGCM outputs and assimilation products, observation-based subduction studies still have been limited to its climatological aspects (Huang and Qiu 1994; Suga et al. 2008) because of the highly restricted spatiotemporal coverage of observational data of the ocean. The situation has recently been changed in large part because of Argo, a global array of over 3500 free-drifting autonomous floats measuring the temperature and salinity of the upper 2000 m of the ocean (Roemmich et al. 2009; Freeland et al. 2010).
The purpose of this study is to demonstrate spatiotemporal variations of the subduction rate based mainly on the Argo observational data. We estimate the subduction rate over the North Pacific based on one of the monthly gridded Argo data products and describe spatiotemporal variation in the subduction rate from 2005 to 2012. Linkages between the subduction rate variability and the pycnocline variability, such as variability in the mode waters and the pycnocline volume, are also discussed.
2. Data and methods
The temperature and salinity data from the grid point value of the monthly objective analysis using the Argo data (MOAA GPV; Hosoda et al. 2008) from January 2004 to February 2014 are used. MOAA GPV is a global monthly dataset of temperature and salinity created using a 2D optimal interpolation method on pressure surfaces. The data are provided at standard depth levels with a horizontal resolution of 1° × 1°. Data from January 2001 to the present are available; however, the data before 2003 are not used in this study since the spatial coverage of Argo profiles in the North Pacific before 2003 is not sufficient.
Figure 1 shows the sparseness of observation data at the 10-dbar depth in February and March each year used for the optimal interpolation. The sparseness is defined as the average distance to N points (N = 4 in Fig. 1) nearest the observation data (Robbins 2014). The sparseness is very large over a quite large area centered at the central to eastern part of the basin especially from 2001 to 2003. Although some patchy data-sparse regions still exist in the central part of the basin even in recent years, the sparseness is generally less than a few hundred kilometers after 2005 in the subtropical gyre. The MOAA GPV dataset produced by the optimal interpolation method reflects the variations in the observation data distribution (i.e., sparseness). The interpolated fields in the large sparseness regions are not so much different from the climatology used as a first guess. The data distribution thus affects the winter MLD and the subduction rate. For example, when a patch of deep MLD maximum in the CMW formation region deeper than the climatology exists but is not observed by any floats in a specific winter, the MOAA GPV will underestimate the winter MLD in the CMW formation region for the year and will cause substantial error in the CMW subduction rate estimate by weakening the horizontal MLD gradient and/or by reducing MLD change along the Lagrangian trajectory in two successive winters. The sparseness map suggests the subtropical gyre region where deep winter MLDs develop is well observed after 2005.
Another point to be noted is its large decorrelation scale used to construct the MOAA GPV dataset. The decorrelation radii of White (1995) were used, whose ranges at the surface in the North Pacific are 10.4°–24.0° in zonal and 6.6°–11.5° in meridional directions. The mode waters in the North Pacific are formed near oceanic fronts and have narrow and patchy structures in their formation regions (e.g., Rainville et al. 2007; Oka et al. 2012). These structures with relatively small horizontal scales including mesoscale or with sharp horizontal gradients are presumably almost smoothed out in the MOAA GPV fields owing to its interpolation procedures. In general, the smoothing makes the MLD gradient less sharp, which in turn reduces the subduction rate.
Temperature and salinity profiles at each grid point are interpolated onto a 10-dbar vertical grid using the Akima spline (Akima 1970). Isopycnal gridded data of pressure, potential temperature θ, specific volume anomaly δ, and dynamic height relative to 2000 dbar with 1° × 1° horizontal and 0.1 potential density σθ vertical grids are then generated using HydroBase 2 software (Curry 2001). Specific volume anomaly and dynamic height are used to calculate the pressure anomaly streamfunction, defined by Zhang and Hogg (1992) on an isopycnal σ:
where p is the pressure, and p′ is the pressure anomaly on the surface σ. The pressure anomaly is defined as , where is the lateral mean pressure on the surface σ over a domain of interest. In this study, the domain (the North Pacific) is defined as the box domain bounded by the 0° and 60°N parallels and the 120°E and 80°W meridians (Fig. 2). The marginal seas (the Okhotsk Sea and Japan Sea) and shallow coastal seas are not included in the analysis. The pressure anomaly streamfunction is converted to isopycnal velocities assuming a zero velocity surface at 2000 m. Trajectories on isopycnal surfaces of particles released at the base of the late winter mixed layer were computed from the obtained velocity fields.
where wEk is the Ekman pumping velocity, υm is the meridional velocity in the mixed layer, hm,0 and hm,1 denote the MLD in the first and second winters at its released and destination locations, respectively, along a Lagrangian trajectory, and the overbar indicates an average over a 1-yr (forward) Lagrangian trajectory calculated based on the velocity fields mentioned above. The second term on the right hand side (rhs) of Eq. (2) starting with β indicates the vertical velocity reduction due to the meridional velocity υ in the mixed layer. The first two terms on the rhs of the Eq. (2) indicate vertical pumping at the base of the winter mixed layer, and the third term indicates lateral induction due to the slope of the mixed layer base. Positive values mean transport downward. Grid points with negative values are regarded as grid points where no subduction occurs.
The transfer of water in the opposite direction, from the permanent pycnocline to the surface mixed layer, is sometimes called obduction (e.g., Qiu and Huang 1995). The annual obduction rate Oann is calculated in almost the same way as Sann:
The overbar indicates an average over a 1-yr Lagrangian trajectory as Eq. (2), but the trajectory here is calculated backward. In this study, the subduction rate Sann is taken to be positive and the obduction rate Oann is taken to be negative unless otherwise noted.
The annual-mean Ekman pumping velocities in Eqs. (2) and (3) were calculated using the wind stress data from Japanese 25-yr Reanalysis Project (JRA-25)/Japan Meteorological Agency Climate Data Assimilation (JCDAS; Onogi et al. 2007). Monthly data from January 2004 to March 2013 are used. The original data with T106 Gaussian grid were linearly interpolated onto 1° × 1° grid to fit in with the ocean data.
The depth at which the potential density increases by 0.125 kg m−3 from that at 10 dbar is defined as MLD (Levitus 1982). Mixed layer density, temperature, and salinity were defined as the averaged values of those within the mixed layer. The late winter MLD values used in Eqs. (2) and (3) were calculated as the average of those in February and March, as is the case of Suga et al. (2008).
3. Spatiotemporal variation of the subduction rate in the North Pacific
First, spatiotemporal variation of the annual subduction rate from 2005 to 2012 is examined. Remarkable year-to-year variations in the subduction rate distribution are apparent (Fig. 3). It is also worth noting that the subduction rate distribution in any individual year is not merely a modest perturbation to the climatological distribution (cf. Fig. 3d and/or Fig. 4 of Suga et al. 2008). Each map shows subduction maxima with grid values of order 150 m yr−1 as the climatological distribution, with the locations where the maxima appear being largely different from year to year.
The yearly spatial distribution patterns of the subduction rate correspond well with those for lateral induction (Fig. 4), a component consisting of the subduction rate together with the vertical pumping component (Fig. 5). The lateral induction distribution for each year shows large spatial variability with some local maxima. The locations of the local lateral induction maxima largely coincide with those of the subduction rate maxima. Furthermore, they also correspond well with the locations of the deep late winter MLD (Fig. 4). Spatial distributions of the vertical pumping, on the other hand, reveal a homogeneous distribution throughout the subtropical gyre with relatively larger values in the eastern and southern part of the gyre and do not resemble the net subduction rate distributions. This suggests that the spatiotemporal variations of the deep late winter MLD largely determines the subduction rate variations through the lateral induction variations as suggested in previous studies (e.g., Liu and Huang 2012). Further examination on the dominance of the MLD variations in determining the subduction rate variations is performed later in this section.
The obduction rate distribution shows large year-to-year variation as does the subduction rate, but the spatial extent where the obduction occurs is relatively confined near the deep MLD region (Figs. 6, 7). Note that the obduction rate is taken to be positive in this paragraph and for Figs. 6, 7, and 8. The obduction rate shows some local maxima exceeding 100 m yr−1 with values only occasionally and very locally exceeding 150 m yr−1. These spatial distribution patterns largely reflect those of both the lateral induction and the late winter MLD (Fig. 7), as is the case of the subduction rate. Contribution from the vertical pumping is small (Fig. 8). Note that the Kuroshio Extension (KE) region around 35°N west of the date line is represented as blank. The obduction rates at these grid points are unable to be calculated because the 1-yr backward trajectories of the released particles are unavailable because of the intersection of the contours of the streamfunction and the coastline of Japan (i.e., the lateral boundary). These complications in computing the obduction rate in the western boundary regions and its extension region may cause great uncertainty in the obduction rate so that we mainly focus on the subduction rate in the remaining part of this paper.
Next, we examined temporal variations of the total subduction rate over the subtropical North Pacific. The total subduction rate is calculated by integrating the grid values of the subduction and obduction rates over the subtropical North Pacific (defined by 15°–45°N and 135°E–110°W as indicated by thick black lines in Fig. 2) for the potential density range from 22.0 to 27.0σθ. The time mean of the total subduction rate averaged over the analysis period is 34.7 Sv. This is quite similar to previous studies using climatologies [33.1 Sv in Huang and Qiu (1994); 33.6 or 37.8 Sv in Suga et al. (2008)], but large year-to-year variations are found (Fig. 9). The total subduction rate varies from 25 to 50 Sv with the minimum in 2009 and the maximum in 2010 (Fig. 9a). The basin-integrated obduction rate, on the other hand, exhibits very small temporal variations ranging from 5 to 10 Sv. As a result of the variations in these two components, the temporal change of the rate of net water transfer across the base of the mixed layer (calculated as the sum of the annual subduction and obduction rates), which varies from 15 to 40 Sv, is largely dominated by that of the subduction rate.
A time series of volume change rate of the permanent pycnocline of the subtropical North Pacific is also shown in Fig. 9a. Here, the volume V of the permanent pycnocline is defined as the volume of water below the late winter MLD and above the 27.0σθ isopycnal, and the volume change rate ΔV for a certain year is calculated as the difference between the pycnocline volumes of that year and the next year (positive if the volume increases with time) divided by the seconds in a year to give units of Sverdrups. It is clear that the year-to-year variations of the sum of the subduction and obduction rates are largely consistent with that of the pycnocline water volume change rate. This suggests that the water volume of the permanent pycnocline is largely controlled by the water flux across the base of the mixed layer and hence the subduction rate. The calculation of the obduction rate faces a number of complications and obstacles, as was discussed previously. However, the very consistent variability between the sum of the subduction and obduction rates and the pycnocline volume change rate implies that the temporal variability of the obduction rate is indeed small, as shown in Fig. 9a.
To identify the dominant cause of the temporal variation of the subduction and obduction rates, relative contributions of the lateral induction and the vertical pumping terms are investigated (Fig. 9b). Lateral induction and vertical pumping terms for the subduction rate integrated over the density range from 22.0 to 27.0σθ reveal amplitudes that are quite comparable in their relative contributions to the subduction rate through the analysis period. When averaged over the analysis period, the lateral induction and vertical pumping account for 16.2 and 18.5 Sv, respectively. However, the year-to-year variation of the subduction rate is predominantly determined by that of the lateral induction term, as found by previous studies (e.g., Liu and Huang 2012). As for the obduction rate, the lateral induction term dominates the vertical pumping term throughout the analysis period, and the year-to-year variability of the obduction rate is largely determined by that of the lateral induction term as is the case of the subduction rate.
According to Eq. (2) of the subduction calculations, variations of the lateral induction term depend on those of both the late winter MLD and the velocity field on isopycnals along which the trajectory of water parcels is to be tracked. We next examine the relative importance of the late winter MLD variations and the velocity (pressure anomaly streamfunction) variations in determining the lateral induction variations by conducting a set of experiments, in a similar way with Liu and Huang (2012). We calculated two types of yearly annual subduction rates by replacing one of the late winter MLD fields or the velocity fields to its climatology for the analysis period. Hereafter, the experiment using climatological MLD is referred to as the MLDclim experiment and the other using climatological velocity is referred to as the Vclim experiment. By comparing temporal variation of the resultant subduction rates, the term responsible for temporal variations of lateral induction can be identified. The time series of integrated lateral induction over the density range of 22.0–27.0σθ from the two experiments reveal a clear contrast; the lateral induction term of the Vclim experiment exhibits large year-to-year variations with similar spatial distribution shown in Fig. 3, while on the other hand, that of the MLDclim experiment indicates little year-to-year variation (figure not shown). The results from these experiments clearly show the dominant role of the late winter MLD variation in determining the spatiotemporal variation of the subduction rate through the lateral induction term.
Another interesting noteworthy point is the highly correlated variability of the subduction rate with the PDO index (Fig. 9c). A correlation between the North Pacific subduction rate and the PDO has been already reported by Qu and Chen (2009) through use of the OFES outputs, but this is the first study to confirm this with observationally based analyses. We will discuss this point in more detail in a later section.
4. Subduction rate for density classes and its relation with mode waters
Variability of the annual subduction rate is examined more closely by dividing the grid subduction values into 0.1σθ density bins for the density range from 22.0 to 27.0σθ. The subduction rate within the subtropical North Pacific bounded by 15° and 45°N parallels and 135°E and 110°W meridians is considered. It is clear that not only the total subduction volume but also its partition into density classes vary considerably from year to year (Fig. 10). As is the case of the spatial distribution map, the individual years can reveal behavior that is distinct from the climatological density class distribution (cf. Fig. 5b of Suga et al. 2008). In the climatological view of subduction, two peaks of subduction appear around 25.3 and 26.2σθ, from both observational (Suga et al. 2008) and model studies (Qu and Chen 2009), which correspond to W/ESTMW and CMW, respectively. The lighter subduction peak appears almost every year during the analysis period, but the denser peak does not. The most notable change occurs at 25.5–26.6σθ (which correspond to CMW density range) where subduction almost completely disappears in 2009, while it is vigorously active in 2005 and 2010. Modes of the subduction rate appear around 25.0σθ (which corresponds to the STMW density) for every year with a maximum value of 3.65 Sv at 25.05σθ in 2009, but the mode density class and the magnitude of the subduction rate at that density class changes considerably from year to year. In 2008 and 2010, a mode of subduction appears at 24.75σθ and 24.45σθ bins with 2.50 and 2.51 Sv, which are lighter than usual 25.05–25.25σθ bins.
To explain the cause of these considerable year-to-year changes in partitioning of the subduction rate into density classes, we further investigate the subduction rate by dividing it by geographical locations (western and eastern regions, as shown by black dashed lines and white lines in Fig. 2, respectively), paying special attention to its relationship with the North Pacific mode waters. The eastern region includes the ESTMW formation region, and the western region includes the WSTMW, CMW of both lighter and denser variety, and TRMW formation regions. The contributions from the two STMWs (WSTMW and ESTMW), which partly share the same density classes, can be separated by dividing the subtropical North Pacific into the two domains.
A histogram of the subduction rate in the eastern region (Fig. 11a) indicates that most subduction occurs at density ranges between 24.05 and 25.05σθ. The mode density class increases from 24.35σθ in 2006 to 25.05σθ in 2009 and then dropped to become 24.35σθ in 2010. Neither 2005 nor 2011 reveal a clear peak density class, and the subduction rate in none of the density classes surpass 1 Sv. Contribution from the lateral induction component for these 2 yr is very small, suggesting late winter MLDs in these years are shallower than the other years. This is easily confirmed from the late winter MLD map (Fig. 4). Maximum MLD in the eastern region over 2006 to 2010 exceeds 150 m; on the other hand, it is less than 150 m in 2005 and 2011, apparently shallower than the other years. For these 2 yr, the vertical pumping component dominates the subduction rates. Since the vertical pumping distribution (Fig. 5) corresponds not to the surface MLD pattern but rather to a larger-scale atmospheric wind pattern, the resulting density class distribution of the vertical pumping (and hence the subduction) shows relatively homogeneous partitioning among density classes, rather than any localized subduction maximum at a specific density class. The MLD in the eastern region in 2009 is the deepest among others, when the mode density class is densest (Fig. 11a) and the integrated subduction rate within the 24.55–25.05σθ range in the eastern region is largest over the analysis period (Fig. 12a).
In the western region (Fig. 11b), the shape of the histogram varies considerably from year to year as that in the whole subtropical North Pacific, with a relatively clear peak manifesting itself at approximately 25.05–25.25σθ almost every year over the analysis period. Since this density class is the same as that of WSTMW, it is suggested that WSTMW is formed every year; however, the formation volume is considerably different from year to year. Subduction volume in a mode density class (within WSTMW mode water density range) is large in 2005, with this gradually decreasing until 2009 and then increasing to reach the largest value of 2.2 Sv at 25.05σθ in 2011. There appears to be another peak at around 26.05–26.35σθ; however, this peak cannot be found in either 2006 or 2009. This density range corresponds largely to that of the L-CMW, suggesting occasional nonoccurrence of the L-CMW formation. Instead, a relatively weak subduction peak appears at around 26.5–26.6σθ in 2006 and 2009. Subduction rates in those densest density classes, corresponding to TRMW/D-CMW density, are largest in 2006.
Figure 12 shows the time series of subduction rates integrated within each mode water density range over each region. The subduction rate for nominal ESTMW (integrated in the east within 24.55–25.05σθ range) varies from about 3.5 Sv in 2005 to 8.8 Sv in 2009 (Fig. 12a). Nominal WSTMW subduction (integrated in the west within 25.05–25.55σθ range) exhibits a minimum value of 1.8 Sv in 2009 and a maximum of 9.9 Sv in 2011 (Fig. 12b). Variation of the nominal CMW subduction rate (integrated in the west within 25.85–26.35σθ range) is larger, varying from 0.5 Sv in 2009 to 11.5 Sv in 2010 (Fig. 12c). The ranges of variability for WSTMW and CMW subduction are larger than that for ESTMW, reflecting larger late winter MLD variability in the WSTMW and CMW formation region in the western region (Fig. 4). The volume change rate of the main pycnocline is also shown in Fig. 12, which is calculated in the same region and in the same density range. The temporal variability of the subduction rate for each mode water corresponds relatively well to that of the pycnocline volume change rate, suggesting small variations in exchange through lateral boundaries and in diapycnal transports of water through the upper and lower isopycnal surfaces.
The time series of the subduction rate for WSTMW and CMW and for ESTMW show anticorrelated decadal variability as suggested by Chen et al. (2010) using OFES output. During the warm phase of PDO, the sea surface temperature in the western (eastern) subtropical North Pacific is cooler (warmer) and hence the mixed layer is deeper (shallower) than usual. The deepening (shoaling) of the winter MLD results in increased (decreased) lateral induction and hence increased subduction rate in the west (east) and vice versa for the cool phase of PDO. If this mechanism works, an increased subduction rate is related not only to the deeper winter MLD but also to the cooler (denser) winter mixed layer. In fact, the mode density range in the eastern region becomes denser (Fig. 11a) as the subduction rate increases from 2005 to 2009 (Fig. 12a) and then become lighter as the subduction rate decreases from 2009 to 2011. Although it is somewhat obscured, a similar relationship between the subduction rate and the mode density class can also be found in the western region within the WSTMW density range (Figs. 11b and 12b). That is, in 2009 when the subduction rate is the smallest over the analysis period, the mode density class of the WSTMW is the lightest. Sugimoto et al. (2013) reported marked freshening of the WSTMW in 2009 and 2010, through excess rainfall in the 2008 warm season. It is obvious that not only temperature but also salinity contribute to changes in density. When considering in detail the relationship among the subduction rate, the water properties of subducted water, and the PDO, it is necessary to take into account each process contribution to the subduction process such as eddy activity, heat and freshwater fluxes, and so forth.
The delayed response of physical conditions in the upstream KE region (WSTMW and western part of L-CMW formation region) to the PDO through the decadal variations of the KE system (Qiu and Chen 2006) is one of the important processes to be considered; the WSTMW formation is expected to have a lag of a few years from the PDO. Oka et al. (2012) analyzed Argo profiling data taking the delayed response into account and suggested that the WSTMW formation in winter just south of the upstream portion of the KE at 140°–152°E was active (reduced) during the stable (unstable) period of the upstream KE in 2002–05 and 2010/11 (2006–09), and this decadal WSTMW formation variability is anticorrelated with that of the L-CMW just north of the KE. It is interesting that the decadal variability of the WSTMW subduction (Figs. 11b, 12b) is apparently in phase with the variability of WSTMW formation suggested by Oka et al. (2012). This correspondence suggests that the WSTMW variability discussed above is actually a delayed response to the PDO. Chen et al. (2010) showed anticorrelated variability in the subduction rate in the western (including WSTMW and CMW) and eastern North Pacific and suggested surface wind and heat flux are the main controlling factors for the subduction rate of the North Pacific. Their composite analysis between strong and weak western subduction years (Fig. 3 of Chen et al. 2010) clearly shows differences in wind and heat flux in the central North Pacific are large; however, the differences are not so large over the WSTMW formation region in the western North Pacific. This also supports the idea that the oceanic forcing is dominant for the WSTMW subduction variability; negative sea surface height anomalies generated in the eastern North Pacific propagate westward, destabilize KE, and reduce WSTMW subduction with a lag of a few years.
On the other hand, the decadal variability of the CMW subduction (Figs. 11b, 12c) is out of phase to that suggested by Oka et al. (2012). It is obvious from the spatial distribution map of the subduction rate (Fig. 3) that there exists little subduction just north of the KE at 140°–165°E. This region has a strong meridional property gradient associated with KEF and KBF. In addition, it is suggested that mesoscale eddies play two important roles for the L-CMW formation variability in this region: the first is that nonlinear anticyclonic eddies actively detached from KE during the unstable period of the upstream KE bring low potential vorticity water from the southern recirculation gyre into the L-CMW formation region, generate a precondition favorable for deep mixed layer development by reducing the background stratification in the L-CMW formation region, and the second is that these eddies themselves serve as important formation sites of L-CMW (Oka et al. 2012). Water properties and MLD in these frontal regions with high eddy activities would be severely modified artificially because of the large decorrelation scale of MOAA GPV. Also, the late winter MLD in this region is generally shallower than that in the downstream region (Fig. 4). It is implied that L-CMW formed in the upstream KE region is reentrained into the surface mixed layer in the downstream region in the subsequent winter and does not contribute to actual subduction. The temporal variations in the L-CMW density range in Fig. 12c thus mostly reflect MLD and subduction variability in the central part of the basin. As mentioned above, Chen et al. (2010) suggested the MLD and subduction variability in the central North Pacific is largely controlled by the atmospheric forcing anomalies associated with the PDO. The stronger westerlies in the L-CMW formation region during the warm phase of PDO makes the heat loss of the surface ocean larger, which results in deeper winter MLD development in the region during the warm phase of PDO.
5. Estimation of error in the subduction rate calculation induced by choice of late winter month
In this study, late winter MLDs are determined by using both February and March data. However, since the MLD sometimes differs substantially between February and March (with the differences as large as 180 m) largely due to interannual variability in the surface buoyancy forcing and partly due to observational density, the choice of the time interval over which one calculates the average depth of the late winter mixed layer could have a significant effect on the resulting calculated subduction rate. Similarly, the 1-yr average indicated by the overbar in Eq. (2) [Eq. (3)] is calculated for the 12 months from March (February) of a certain year to February (March) of the next (previous) year. The selection of month from which the 1-yr average starts might also have nonnegligible consequences for the subduction rate as well.
To estimate the amplitude of the errors mentioned above, we conducted nine types of the subduction calculations with different combinations of three types of late winter mixed layer data (derived from February data only, March data only, and average of February and March) and three types of 1-yr-averaged hydrographic fields (with 1-yr average starting from January, February, and March). The time series of integrated subduction rates over the subtropical North Pacific resulting from the nine calculations are shown in Fig. 13. Using the February and March averaged mixed layer data with different 1-yr-averaged hydrographic fields results in approximately a 2.2% error in the subduction rate integrated over the subtropical North Pacific on average over the period from 2005 to 2012. The error itself varies substantially from year to year, reaching as large as 4.9%. Using the 1-yr-averaged hydrographic data starting in March, with different mixed layer data, results in an error of approximately 5.0%. This error also largely varies from year to year and is as large as 10.6%. The selection of the time period over which the late winter mixed layer properties are to be calculated causes a larger error in the resultant subduction and obduction rates than the selection of the month from which the 1-yr average in Eqs. (2) and (3) start. As mentioned earlier, the observational data distribution would also cause significant error in representing the winter MLD and the subduction rate. For a proper estimate of the subduction rate, two successive winter MLDs along the Lagrangian trajectory should be property reproduced. According to the sparseness map (Fig. 1), the winter MLD in the subtropical gyre is well observed after 2005 by using both February and March data. We speculate that the error caused by any limitation in the data distribution might be in similar order with the error estimated by using different winter MLDs calculated only using February or March data. The errors are not negligible when considering the small temporal undulation of the subduction rate; however, the discussion about the interannual to longer time-scale variability, with magnitude of as large as 50% of the total, demonstrated in this study is still valid.
Through analysis of a monthly gridded Argo data product, the spatiotemporal variability of the subduction rate in the North Pacific during 2005–12 has been evaluated based mainly on observation data in this study. Substantial interannual variability has been found both in the spatial distribution of subduction as well as in the area-integrated subduction rate. The total subduction rate integrated over the subtropical North Pacific varies from 25 to 50 Sv, with a time-mean value of 34.7 Sv over the analysis period. The spatiotemporal variability of the subduction rate is largely controlled by that of the late winter MLD through the lateral induction term of the subduction rate, though the relative contribution of the lateral induction and the vertical pumping terms are almost comparable on average. The temporal variability of the integrated subduction rate over the subtropical North Pacific is well correlated with the volume change rate of the permanent pycnocline water. The integrated subduction rate over the CMW and ESTMW formation regions reveal anticorrelated temporal variability, which is simultaneously correlated with the PDO. The WSTMW subduction rate variability is presumably correlated with the PDO with a lag of a few years. Partitioning of the subduction rate into 0.1 kg m−3 density classes reveals that the increase (decrease) in the subduction rate of both the WSTMW and ESTMW area accompanies densification (lightening) of the mode density class. The subduction rate in the CMW density classes also undergoes large interannual variation; it is robust and large in 2005 and 2010 but almost disappears in 2009.
The smoothing effect of the interpolation method used to produce the gridded Argo product on the estimation of the subduction rate is not trivial. In general, smoothing procedures tend to make the horizontal MLD patterns less patchy and its gradient less sharp and thus tend to underestimate the subduction rate in the real ocean, which may be more variable both in time and space associated with meso- to smaller-scale phenomena. The calculation error estimated in this study thus may be an underestimate; however, the general consistency with previous modeling and observation results lends credibility to our analysis. As it is clearly shown that the subduction rate variability has great influence on the water volume changes of the main pycnocline, one can easily expect that the subduction rate variability will also have significant influence not only on the variability of the subsurface physical properties such as potential vorticity, but also on that of the biogeochemical fields such as carbon, oxygen, and nutrients. Combined use of the observational data with ocean GCM and/or Earth system models will help us better understand the physical and biogeochemical variability of the upper ocean associated with variability in the physical subduction processes.
The authors appreciate the valuable comments offered by members of the Physical Oceanography Group at Tohoku University, Keith Rodgers, and two anonymous reviewers. The first author (KT) was supported by the Global Center-Of-Excellence (GCOE) Program, “Global Education and Research Center for Earth and Planetary Dynamics” at Tohoku University. The second author (AI) was supported by the Sasakawa Scientific Research Grant from Japan Science Society. The third author (TS) was supported in part by the funds from the Japanese Society for Promotion of Science [Grants-in-Aid for Scientific Research (B) 21340133 and 25287118] and from the Ministry of Education, Culture, Sports, Science and Technology [Grants-in-Aid for Scientific Research on Innovative Areas 22106007 (“A ‘hot spot’ in the climate system: Extra-tropical air–sea interaction under the East Asian monsoon system”)]. This work was also supported in part by NSF Award OCE-1155983.
Current affiliation: Oceanography and Geochemistry Research Department, Meteorological Research Institute, Tsukuba, Japan.