Ocean Acidification State in the Highly Eutrophic Tokyo Bay, Japan: Controls on Seasonal and Interannual Variability

Seasonal and interannual variabilities in the partial pressure of CO2 (pCO2), pH, and calcium carbonate saturation state (Ω) were investigated in the highly eutrophicated Tokyo Bay, Japan, based on monthly observations that were conducted from 2011 to 2017. There were large variabilities in these parameters for surface and bottom waters due to photosynthesis and respiration, respectively. Warming/cooling and freshwater input also altered the surface Ω. During the observation period, calcium carbonate undersaturation was observed twice in the anoxic bottom waters in summer. The data indicate that anaerobic remineralization under anoxic conditions lowers the Ω, causing undersaturation. These findings suggest that de-eutrophication can decelerate ocean acidification in the bottom waters of Tokyo Bay. However, if atmospheric CO2 exceeds 650 ppm, aragonite undersaturation will be a common feature in the summer bottom water, even if hypoxia/anoxia are alleviated by de-eutrophication.


INTRODUCTION
Ocean acidification (OA) is defined as the progressive reduction in the pH of ocean over an extended period of time, typically decades or longer (Field et al., 2011). Changes of seawater carbonate chemistry associated to OA, such as an increase in partial pressure of seawater CO 2 (pCO 2 ) concentration and a decrease in calcium carbonate (CaCO 3 ) saturation state ( ) are also included in the discussion of the effects of OA (Royal Society, 2005;IPCC, 2014). The increased CO 2 concentration promotes the photosynthesis by primary producers such as phytoplankton and seagrass (Hinga, 2002;Palacios and Zimmerman, 2007;Hendriks et al., 2010), while the decreased adversely affects the calcification, as well as the growth and survival of calcifying organisms (e.g., Kroeker et al., 2010). The increased CO 2 concentration also affects the growth of otoliths (CaCO 3 structures in ears), metabolic rate and behavior of some fishes (Heuer and Grosell, 2014 for review). The co-occurrence of OA with other environmental changes, such as warming and deoxygenation lead to further biological complications (Kroeker et al., 2013;Gao et al., 2019;Doney et al., 2020). It has been found that an elevated CO 2 concentration reduces the tolerance of fish and bivalves to low oxygen conditions (Hancock and Place, 2016;Miller et al., 2016;Gobler and Baumann, 2016). In addition, the synergistic effects of OA and the warming of seawater negatively impact the calcification, reproduction, and survival of marine organisms (Harvey et al., 2013).
The current OA, caused primarily by the uptake of anthropogenic CO 2 from the atmosphere, is occurring at a rate higher than that observed at any other time in the last 300 million years (Hönisch et al., 2012). The high OA rate, together with other anthropogenic perturbations, pose a threat to marine organisms, and thus to human society (Bednaršek et al., 2016;Hall-Spencer and Harvey, 2019;Stewart-Sinclair et al., 2020). Coastal regions are of particular importance to human society, providing ∼40% of the global ecosystem goods and services (Costanza et al., 1997). However, due to the complex natural processes and anthropogenic effects, it is difficult to understand and predict the progress of OA in coastal regions. For instance, freshwater discharge from land lowers the salinity and alkalinity and thus lowers the buffer capacity of coastal seawater, making it sensitive to anthropogenic CO 2 addition (Salisbury et al., 2008;Hu and Cai, 2013). An excessive loading of man-made nutrients from land, a phenomenon known as eutrophication, causes massive algal production (blooms) and the subsequent decomposition of large amounts of sinking organic matter causes hypoxia (oxygen concentrations ≤ 2 mg L −1 or ≤ 61 µmol kg −1 ) in coastal waters worldwide (Breitburg et al., 2018). Eutrophication can thwart OA in the surface waters by increasing the biological CO 2 uptake; however, it amplifies OA in bottom waters by increasing the amount of CO 2 released from the decomposed organic matter (respiration) and also by reducing the buffer capacity (CO 2 -rich water has lower buffer capacity; Cai et al., 2011;Duarte et al., 2013;Wallace et al., 2014;Zhao et al., 2020). Large seasonal variations in the physical and biogeochemical conditions of shallow coastal regions can also result in large pH and fluctuations, suggesting that the progress of OA causes coastal waters to seasonally exceed the biologically harmful pH and thresholds faster than offshore waters (Pacella et al., 2018;c.f., McNeil and Matear, 2008). Depending on the balance between the OA (caused by anthropogenic CO 2 ) and other processes, the coastal regions can undergo a variety of changes from rapid acidification to basification (Duarte et al., 2013;Ishizu et al., 2019). Furthermore, all of these natural and anthropogenic processes can differ across regions. Therefore, to predict the future environmental conditions, it is essential to quantitatively understand the processes affecting OA in each coastal region (e.g., Cai et al., 2011;Laurent et al., 2017;Feely et al., 2018;Rheuban et al., 2019).
In this study, we investigate the state of OA in Tokyo Bay, a highly eutrophicated bay in Japan. Our previous study was the first report on the of water in Tokyo Bay based on monthly observations at two stations, which were conducted from April 2011 to January 2012 . Of these two stations, larger seasonal variabilities in the physical and biogeochemical parameters were observed at F3, the inner bay station, than at F6, the middle bay station. Herein, the seasonal and interannual variations in the various carbonate chemistry parameters, viz., CO 2 , pH, and are quantitatively examined to understand the processes affecting the of OA in Tokyo Bay, based on the monthly observations made at F3 from April 2011 to December 2017.

STUDY AREA
Tokyo Bay is a semi-enclosed bay with a mean depth of 19 m (Figure 1). The freshwater discharge from rivers (∼400 m 3 s −1 ) and the inflow of offshore waters drives the estuarine circulation in the bay. The residence time of water in the bay is 1-3 months (Okada et al., 2011). The environment of Tokyo Bay has significantly changed over the last century (see Furukawa and Okada, 2006 for review). The nutrient loading from the surrounding metropolitan areas to the bay rapidly increased in the 1950s and the 1960s. The massive algal bloom (commonly referred to as the "red tide"), mainly caused by the dinoflagellates and diatoms, and the resulting bottom water hypoxia/anoxia have been repeatedly observed since the 1960s Kodama and Horiguchi, 2011). Since most river discharge enters from the north-west side, severe eutrophication and prolonged hypoxia/anoxia are observed in the inner part of the bay (e.g., Kodama and Horiguchi, 2011). Hypoxia/anoxia promotes the anaerobic oxidation of organic matter such as denitrification and sulfate reduction (e.g., Krumins et al., 2013). The bottom water with toxic hydrogen sulfide produced by the sulfate-reducing bacteria can rise by wind mixing, killing the fish and shellfish at the northern coast of the bay. This phenomenon is known as "blue tide" because the hydrogen sulfide in upwelled water is oxidized to elemental sulfur, causing the water to appear whitish blue (Takii et al., 2002). In addition, the area of the tidal flats decreased by ∼90% over the past century, thereby decreasing the shellfish habitat and destroying the network for benthic larvae transport. Due to these environmental changes, the total fish caught annually in Tokyo Bay decreased from > 1 × 10 5 tons in the 1950s (∼90% of shellfish) to 4 × 10 4 tons in the 1990s (Furukawa and Okada, 2006). Although the implementation of sewage treatment decreased the nutrient concentration in the bay after the 1980s (Kanda et al., 2008;Kubo et al., 2019), red tide and summer bottom hypoxia are still observed every year. The progress of OA can impose additional stress to the ecosystem of Tokyo Bay, which has been suffering from eutrophication and hypoxia.

MATERIALS AND METHODS
Monthly observations (a total of 77 times) were carried out at F3, a time-series station (35.51 • N, 139.83 • E, bottom depth = 23 m; Figure 1) from April 2011 to December 2017. Station F3 is located in the inner bay region. Seawater samples were collected using Niskin bottles mounted on a conductivitytemperature-depth (CTD)/rosette on the R/V Seiyo-maru of Tokyo University of Marine Science and Technology at target depths of 0 m, 5 m, 10 m, 15 m, and bottom. Actual depths of sampling were 1-3 m for the target depth of 0 m, 2-8 m for 5 m, 10-14 m for 10 m, 14-18 m for 15 m, and 19-24 m for the bottom.
Samples for dissolved inorganic carbon (DIC) and total alkalinity (TA) were preserved using a saturated mercuric chloride solution (0.1% by volume). The DIC was determined using the CO 2 extraction-coulometry method . The TA was measured using a spectrophotometric method (Yao and Byrne, 1998) for samples from April 2011 to January 2012, and by an open titration method  for all other samples. The DIC and TA were calibrated against certified reference materials provided by Dr. Dickson at Scripps Institution of Oceanography and/or KANSO Technos, Japan. The pooled standard deviations obtained from replicate analysis were 2.5 and 2.8 µmol kg −1 for TA and DIC, respectively. The dissolved oxygen (DO) and chlorophyll a (Chl a) concentrations were measured using a RINKO DO sensor (JEFF Advantec Co., Ltd.) and a Seapoint chlorophyll fluorometer (Seapoint sensors, Inc.), respectively. The percentage of oxygen saturation (DO%) was calculated at standard air pressure with solubility of Weiss (1970). We also used nutrient data till October 2015, reported in Kubo et al. (2019).
The DIC, TA, temperature (T), salinity (S), and phosphate and silicate concentrations were used to calculate the in situ partial pressure of CO 2 (pCO 2 ), pH (total scale), and of the collected samples using the CO2SYS program (Lewis and Wallace, 1998), along with the K 1 and K 2 constants reported by Mehrbach et al. (1973), refit by Dickson and Millero (1987), and K SO4 constant of Dickson (1990). The phosphate and silicate concentrations of Kubo et al. (2019), who reported data till October 2015, were used. After October 2015, we used the mean of the monthly phosphate and silicate concentrations calculated from 2011 to 2015 at each depth. For , only the saturation state of aragonite ( ar) is reported. The calcite saturation state is approximately 1.5 times higher than that of ar. Figure 2 shows mean values of physical and biogeochemical parameters for each month and each target depth (0, 5, 10, 15 m, and bottom). In December and January, all observed parameters show relatively uniform values from the surface to the bottom, because of vertical mixing induced by cooling and winds. From January to August, the T increases from 10.7 to 26.8 • C in the surface layer and from 12.5 to 21.2 • C in bottom layer (Figure 2A). The S below 10 m remains largely unchanged throughout the year, ranging between 33 and 34 at the bottom ( Figure 2B). In contrast, there is a large seasonal variability in S at 0 and 5 m, decreases to 26.1 and 29.7, respectively, in September. This indicate that river discharge significantly lowers S in the top 5 m layer at F3. The average T and S for all observations were 17.5 • C and 32.0, respectively.

Seasonal Variability
Chl a concentrations in the top 5 m layer also show large seasonal variability ( Figure 2C). The highest Chl a concentration of 28 µg L −1 was observed in the surface water in May. Interestingly, the Chl a concentration increased from December to January at all depths ( Figure 2C), suggesting the enhancement of primary production during January. In addition, the concurrent decrease in the phosphate and silicate concentrations from December to January (not shown) indicates the occurrence of photosynthesis by diatoms. This finding is consistent with previous studies, which have reported winter photosynthesis by diatoms in Tokyo Bay, stimulated by episodic stratification with favorable light conditions and abundant nutrients (Shibata and Aruga, 1982;Nomura and Yoshida, 1997).
The surface DO was slightly undersaturated in December (94.5 ± 3.9%), but increased to 115.0 ± 15.5% in January due to the production of DO by photosynthesis. Despite continued photosynthesis and warming (lowers gas solubility in seawater), DO% does not largely increase from spring to summer, possibly due to the release of O 2 to the atmosphere by an air-sea gas exchange. From spring to summer, the DO in the bottom water at F3 largely decreases because of decomposition of organic matter, reaching as low as 3.3% (7.4 µmol kg −1 ) in August. Bottom water is in hypoxic condition (≤61 µmol kg −1 ) from July to September. In August, water at 15 m also reaches hypoxic condition. After September, vertical mixing supplies oxygen from the atmosphere to the bottom waters, increasing the DO levels ( Figure 2D).
Carbonate chemistry should be affected by above mentioned seasonal variation in cooling/warming, freshwater mixing, biological activity such as photosynthesis and respiration. Seasonal variations in TA and DIC were similar to those in S, showing low values during summer in the top 5 m layer due to mixing with river water (Figures 2B,E,F). The DIC below 10 m increases from spring to summer due to the accumulation of CO 2 from decomposition of organic matter. The seasonal and interannual variations in the carbonate parameters (calculated from TA and DIC) differ from those in TA or DIC (Figures 2G-I). The seasonal variation in pCO 2 is much larger in the bottom layer than in the surface layer. Surface water tends to be undersaturated with respect to atmospheric CO 2 (∼400 µatm) from January to October, while it is oversaturated with CO 2 from November to December ( Figure 2G). The bottom water also tends to be undersaturated from January to April and the pCO 2 in the bottom water increased as the DO% decreased, reaching 1,134 µatm in August. From September to November, the bottom pCO 2 rapidly decreases, indicating the transfer of CO 2 to the upper waters and the atmosphere by an enhanced vertical mixing in fall; however, CO 2 oversaturation was still observed in December, following which, winter photosynthesis decreases the pCO 2 till the water is undersaturated with respect to the atmosphere in January ( Figure 2G). The seasonal variations in the pH mirrors that of pCO 2 ( Figure 2H). Bottom pH is as low as 7.6 in August.
ar is high in the surface waters during summer ( Figure 2I), in spite of freshwater mixing that lowers ar (Salisbury et al., 2008;Yamamoto-Kawai et al., 2011). This means that effect of biological activity (i.e., uptake of CO 2 by photosynthesis) exceeds the effect of freshening. The bottom ar reaches a minimum value of 1.3 in August.

Interannual Variability
Although 7 years of observations are not sufficient to determine the acidification trends in the highly variable Tokyo Bay, the interannual variabilities in the physical and biogeochemical parameters can be investigated. Interannual variability in T is not so large, compared to seasonal variability ( Figure 3A).  (Figure 2A). In the layer below 10 m depth, S is relatively constant throughout the study period. Accordingly, extreme surface freshening leads to strong stratification at around 5 m.
Though a seasonal maximum Chl a concentration of 28 µg L −1 was found in May in Figure 2C, greater concentrations were observed during other months too ( Figure 3C). According to Nomura (1998), a Chl a concentration > 30 µg L −1 indicates the occurrence of red tide in Tokyo Bay, and such high concentrations were observed every year from 2011 to 2014, with the highest Chl a concentration (92 µg L −1 ) in June 2014. The observed Chl a concentrations were lower than 30 µg L −1 after 2015. The extent of DO oversaturation was extremely high (>150%) in the surface layers in 2012, 2013, and 2014. These observations along with the high Chl a concentrations suggest that the primary productivity was higher in the earlier years of the observation period (Figures 3C,D). On the contrary, a significant undersaturation of DO in surface water (<80%, down to 47%) was observed in July and August 2015, October 2016, and July 2017, indicating the mixing with waters affected by respiration at the bottom of F3 or somewhere in upstream regions. In October 2016, the low DO% in extremely low S water (S = 16.3) suggests that the respiratory processes took place near the river mouth and then advected to F3.
Bottom hypoxia (DO ≤ 61 µmol kg −1 , DO% < 25%) was observed every summer. Anoxic conditions (DO = 0 µmol kg −1 ) were observed in the bottom layer in August 2013, July 2015, and August 2015, and at a depth of 16 m in July 2014 ( Figure 3D). Not anoxia, but very low concentration of DO (<5.0 µmol kg −1 , DO% < 2%) were observed in the bottom water in all other years. Although hypoxia and anoxia are induced by the stratification and decomposition of organic matter produced by primary production, the occurrence of anoxia at Frontiers in Marine Science | www.frontiersin.org the bottom does not seem to be directly related to the strength of stratification (stronger in 2011 and 2016 by freshening; Figure 3B) or surface productivity (higher in 2011-2014 than 2015-2017; Figure 3C). This is likely because bottom anoxia was not a special event for 2013, 2014, and 2015, but almost occurs every year, and our observations happened to be made during a short anoxic period in 2013, 2014, and 2015. This is supported by the data of public water quality measurement near our study site (at 9 stations in the area of 35. 46-35.53 • N, 139.74-139.85 • E; data can be downloaded from websites of Tokyo Metropolitan Government Bureau of Environment and Kanagawa prefecture Environment and Agriculture Bureau). The minimum DO values at or below detection limit [<0.5 mg L −1 (Tokyo) or 0.1 mg L −1 (Kanagawa)] were reported every summer at 1-4 of 9 stations, suggesting possible occurrence of localized anoxia. However, such low DO condition was not observed for two consecutive months at each station. These results suggest that anoxia in this region is not uncommon, but is a localized, short-term event. This is consistent with previous studies reporting that bottom DO in Tokyo Bay varies in a short period of time due to the wind-driven intrusion of offshore water (Fujiwara and Yamada, 2002;Sato et al., 2012). The intrusion of offshore water supplies oxygen and also moves hypoxic/anoxic water toward the coast and to the subsurface layer (Fujiwara and Yamada, 2002). This explains the observation of anoxic water at a depth of 16 m, not at the bottom in July 2014.
Interannual variations in TA and DIC are similar to that of S in the upper layer (Figures 3E, 4F). Especially low TA and DIC were observed in September 2011, April 2014, June 2014, and August 2016 when extreme freshening was observed. Although extreme freshening was also observed in October 2016, TA and DIC were relatively high. In the deeper layer, T shows relatively small interannual variability. The bottom DIC in summer is high in 2013, 2014, 2015, and 2017. Lower pCO 2 , higher pH and ar were observed in the summer surface waters during the first half of the observation period (Figures 3G-I). Extremely high pCO 2 (>1,500 µatm), low pH (<7.5) and low ar (<1) were observed in the bottom layer in August 2013 and at a depth of 16 m in July 2014, when anoxic condition was also observed. This is the first record of aragonite undersaturation in Tokyo Bay. Although bottom anoxia was observed in July and August 2015, the carbonate parameters did not show extreme values and ar was higher than 1 (1.2 and 1.4, respectively). In November and December of 2014, high pCO 2 , low pH, and low ar were observed in the whole water column compared to other years. This should reflect the accumulation of large amounts of CO 2 during the preceding summer in the layer below 10 m (Figure 3G), indicating that the properties of summer bottom water later affect the properties in winter. Another process causing interannual variation in the winter carbonate properties is winter photosynthesis (Shibata and Aruga, 1982;Nomura and Yoshida, 1997): relatively low pCO 2 , high pH and high ar were observed in the whole water column in January 2012 and February 2013 when Chl a was relatively high.

Processes Affecting Carbonate Chemistry in Tokyo Bay
Physical Processes Occurring With Air-Sea CO 2 Equilibrium In this section, the processes responsible for the seasonal and interannual variations in the carbonate parameters in Tokyo Bay are quantitatively discussed. As shown in Figure 4, both DIC and TA are linearly correlated with S, indicating the significant influence of mixing of offshore water (high in S, DIC, and TA) with freshwater (low in S, DIC, and TA). The TA and DIC values were higher than expected from S in October 2016 (circled data in Figure 4) when freshening was observed from the surface (S = 16.3) to a depth of 16 m (S = 28.2 at 16 m) in the water column. Although the reason is not completely understood, the low DO% (68% at the surface) suggests that this anomaly arises from the input of local freshwater input with especially high TA and DIC and local respiratory processes in the regions near the coast or the river mouth. Except for the data recorded in October 2016, TA correlated well with S and the regression line of TA is 37.61 (± 0.73) × S + 998.4 (± 23.1) (black line in Figure 4A, r = 0.94, n = 368), which is in good agreement with that obtained by Taguchi et al. (2009) for surface waters in 2008 (TA = 36.77 × S + 1,006, r = 0.99, n = 59), following a mixing line between offshore water [North Pacific Sub-Tropical Mode Water (NPSTW), S = 34.82, TA = 2293 µmol kg −1 ] and river water (S = 0, TA ∼1,000 µmol kg −1 ). This indicates that the waters in Tokyo Bay can be expressed as a mixture of these two components. Deviations from the regression line in the lowsalinity waters may reflect the variations in TA in the river water and the contribution of precipitation (TA = 0). The correlation between DIC and S, except for data from October 2016, is expressed by DIC = 53.7 (± 1.7) × S + 261.5 (± 54.8) (black line in Figure 4B, r = 0.85, n = 367). Compared to the TA-S plot, there were more data scattering and deviations from the mixing line between the DIC values of offshore waters (DIC = 1,960-1,900 µmol kg −1 at S = 34.82 from Ishii et al., 2011) and river waters (1,047 ± 164 µmol kg −1 in the Tamagawa River and 1,499 ± 222 µmol kg −1 in the Arakawa River from May 2011 to April 2012, Kubo, 2015; ∼1,800 µmol kg −1 in the Tamagawa River in April 1993, Ogawa and Ogura, 1997) in the DIC-S plot. Water samples with high S and low DO% showed an elevated DIC, while surface waters with low S and high DO% showed a lower DIC than that expected from the mixing line, indicating the significant influence of photosynthesis and respiration on the DIC concentrations. The DIC should also have been affected by warming/cooling and air-sea gas exchanges. To separately and quantitatively assess the effects of the physical processes on the carbonate parameters, we have assumed that the TA and DIC values were affected only by the physical processes, i.e., mixing with freshwater and warming/cooling. Under this assumption, the TA (TA phys ) can be predicted from the S values assuming two end-member mixing of offshore water with river water (TA phys = 36.76 × S + 1,013). For the offshore waters, the S and TA values of Taguchi et al. (2009) were used (34.82 and 2293 µmol kg −1 , respectively). This TA value was in agreement with that of Ishii et al. (2011), who reported a mean salinitynormalized TA of 2,305 µmol kg −1 at S = 35, corresponding to 2,293 µmol kg −1 at S = 34.82 in the Kuroshio region, southwest of Tokyo Bay. For the river water, a TA of 1,013 µmol kg −1 , the mean of monthly observations of four major rivers (the Edogawa, Arakawa, Nakagawa, and Tamagawa rivers, ± 156 µmol kg −1 , n = 42), was used (reported in the appendix of the 2018 annual report on water quality by the Bureau of Waterworks Tokyo Metropolitan Government) 1 . The DIC can be also predicted from the S and end-member DIC values. However, the mixing of river water with offshore water causes an oversaturation of pCO 2 in the waters compared with the atmospheric pCO 2 . Accordingly, the DIC should decrease from the river mouth to the observation site by the outgassing of CO 2 , as observed in Tokyo Bay  and other coastal regions (e.g., Zhao et al., 2020). Therefore, we used pCO 2 instead of DIC, by assuming that pCO 2 in the physics only case (pCO 2 phys ) is at equilibrium with atmospheric pCO 2 at TA phys and observed T and S conditions. For this calculation, the monthly mean dry air CO 2 concentration (xCO 2 ) observed at Mauna Loa (NOAA/PMEL) was converted to pCO 2 at 100% humidity, standard air pressure, and the observed water T. Then, the DIC, pH and of the waters, under the where the function f refers to the CO2SYS thermodynamic equations, T obs and S obs are observed T and S. DIC phys shows relatively large variability (ranging from 1,481 to 2,071 µmol kg −1 ) under the equilibration with atmospheric pCO 2 , reflecting significant effects of T, S, and TA on DIC ( Figure 5A). Of these, the effect of T can be estimated by comparing DIC phys with that calculated by replacing T obs by mean T of 17.5 • C. It is estimated that summer warming to 27 • C decreases DIC by ∼70 µmol kg −1 and winter cooling to 11 • C increases DIC by ∼50 µmol kg −1 . This indicates that the large variability in DIC phys is mostly caused by freshening ( Figure 5A). Nevertheless, when compared to the observed seasonal and interannual variabilities in the carbonate parameters, those calculated for the physical processes with an air sea CO 2 equilibrium were significantly lower (Figures 2F-I, 3F-I, 5A-D, note that the scales of the y-axis in Figures 5B,C are much smaller than those in Figures 2G,H, 3G,H). Figure 5C shows that the pH phys reflects seasonal variations and long-term trend of pCO 2 phys . In addition, it is apparent that the pH phys is affected by S. However, even under extreme freshening, the pH phys is lowered by only ∼0.03, suggesting that the effects of mixing of seawater with river water and warming/cooling on the pH are small after sufficient air-sea gas exchange has occurred, even in Tokyo Bay with large seasonal variability in T and S. Seasonal variation in the phys values is relatively large (Figure 5D). In general, phys values are high in the summer (∼2.8) and low in winter (∼1.8), and the data are positively correlated with T (r = 0.86), reflecting the effect of T on solubility of aragonite (Mucci, 1983). In addition, surface freshening lowers the phys by ∼1 (for example, the difference between the phys value at 0 m (S = 15.87) and that at 5 m (S = 29.82) in September 2011 is 1.3). Due to the combined effects of T and freshening, the lowest and second lowest phys values ( phys = 1.3 and 1.4) were observed in the surface waters in April 2014 (S = 20.75, T = 15 • C) and October 2016 (S = 16.3 and T = 22 • C), respectively, not in the freshest but warm water in September 2011 (S = 15.9, T = 26 • C, phys = 1.6). In terms of interannual variability, freshwater mixing is the major factor affecting the phys values since T has a small interannual variation (Figure 2A). These findings suggest that warming/cooling induces the seasonal variation in phys , whereas fluctuations in freshwater mixing cause the interannual variations in the surface phys in Tokyo Bay.

Other Processes
The net effect of processes other than physical processes on the carbonate parameters can be estimated by the differences between the observed values and those calculated under the physical processes only condition; The other processes that affect the carbonate parameters include photosynthesis, remineralization, formation and dissolution of CaCO 3 , and air-sea CO 2 disequilibrium. Uncertainties and fluctuations in the end-member TA values used to estimate the carbonate parameters for the physical processes-only case can affect the DIC, pCO 2 , pH, and ar values ( values). The seasonal and interannual variations in the values are similar to those in the observed carbonate parameters (Figures 2F-I, 3F-I, 5E-H), indicating the importance of non-physical processes in Tokyo Bay. The average net effects of these processes in the surface (bottom) waters are as follows: DIC = −76 ± 102 µmol kg −2 (+ 50 ± 73 µmol kg −2 ), pCO 2 = −130 ± 127 µatm (+ 192 ± 283 µatm), pH = + 0.17 ± 0.17 (−0.12 ± 0.17), and ar = + 1.0 ± 0.9 (−0.4 ± 0.7). The relationships between carbonate parameters ( pCO 2 , pH, and ar) and the DO% (Figure 6) indicate that photosynthesis and aerobic respiration are the major processes controlling the values. In waters with high DO oversaturation, the effect of non-physical processes is high: maximum values ( pCO 2 : −345 µatm; pH: + 0.6; ar: + 3.3) were observed in the water with the highest DO% (174%). Relatively low pCO 2 along with high pH and ar values were observed in the surface waters during the first half of the observation period ( Figure 5). These findings are consistent with the higher photosynthesis rate in the earlier years, as indicated by the higher Chl a concentration and DO% (Figures 3C,D).
During summer, high pCO 2 values (>1,000 µatm) were estimated for the anoxic and aragonite undersaturated bottom waters collected on 19 August, 2013 and 04 July, 2014 ( Figure 5F). An increase in the pCO 2 without a corresponding change in the DO level under anoxic conditions ( Figure 6A) suggests the influence of anaerobic processes, such as denitrification and sulfate reduction in the sediment, which are followed by a DIC efflux to the overlying bottom water. Both denitrification and sulfate reduction are known to occur in the sediments of Tokyo Bay (Koike and Hattori, 1978;  Nishio et al., 1982;Takii et al., 2002). Significant denitrification in the bottom waters is also expected to occur in summer (Koike, 1993). In our observations, relationship between dissolved inorganic nitrogen (DIN: nitrate + nitrite + ammonia) and dissolved inorganic phosphorous (phosphate) indicates effects of anaerobic processes. As shown in Figure 6, N/P ratio in low oxygen water is much lower than the Redfield ratio of 16 (Redfield et al., 1963), although surface water in the bay has much larger N/P ratio (Kubo et al., 2019). The low N/P ratio can be caused by denitrification (removal of DIN) and/or efflux of phosphate from anoxic sediments (Golterman, 2001). Low N/P values are found in low DO% waters overall, not only in anoxic waters (Figure 6). This indicates that anaerobic remineralization, not only in anoxic bottom waters but also in anoxic bottom sediments, influences bottom water properties in the region. The DIN deficit (DIN -P × 16) in summer hypoxic/anoxic bottom water (DO < 61 µmol kg −1 at 20 m in July, August, and September) was −24.0 µmol kg −1 (−38.4 ∼ −9.0 µmol kg −1 ). If this deficit represent the removal of DIN by denitrification, DIC produced by denitrification is estimated to be 21.2 µmol kg −1 (8.0∼33.9 µmol kg −1 ; DIN deficit/(104 + 16) × 106; cf. Gruber and Sarmiento, 1997). Without this DIC, carbonate parameters in summer bottom water would change significantly: pCO 2 decreases by 136 µatm while pH and ar increases by 0.06 and 0.18, respectively. However, given that DIN deficit can be caused by efflux of phosphate from sediments (Golterman, 2001) and that N/P ratios in Tokyo Bay without denitrification should be greater than 16 (Kubo et al., 2019), these estimates could be either underestimates or overestimates. Nevertheless, these results demonstrate the importance of anaerobic remineralization processes on carbonate parameters in Tokyo Bay. Figures 3 and 5 show that anaerobic remineralization in the water and sediments increases pCO 2 to as high as > 1,500 µatm and pCO 2 as high as > 1,000 µatm, leading to sporadic aragonite undersaturation in the bottom waters. Notably, aragonite undersaturation and anoxia were observed at a depth of 16 m, and not in the bottom waters (22 m) in July 2014. The higher ar and DO values in the waters above and below 16 m imply that the aragonite undersaturation and anoxia were not developed in situ at station F3, but they initially occurred in the surrounding areas and then advected to station F3. Such advection of low DO water has been reported as a subsurface hypoxic water mass in a previous study focused on Tokyo Bay (Fujiwara and Yamada, 2002). In July and August 2015, when bottom anoxia was also observed, the ar values were greater than 1. It seems that DIC in these waters were relatively high (∼175 µmol kg −2 ) but were not enough to cause aragonite undersaturation. As DIC should reflect total amount of organic matter remineralized through both aerobic and anaerobic oxidation, anoxic water does not necessarily always accompany aragonite undersaturation.
The scatter in Figure 6 reflect the effects of non-physical processes other than photosynthesis and respiration on the carbonate parameters. For example, the waters collected in December tended to have higher pCO 2 and lower pH and ar values than those that were expected from the DO% (Figure 6), possibly due to the occurrence of an air-sea CO 2 disequilibrium after vertical mixing. Vertical mixing during fall/winter brings CO 2 -enriched and DO-depleted bottom water to the surface. An air-sea gas exchange changes the properties to attain a DO% = 100% and = 0. However, as compared to the DO exchange rate, the CO 2 exchange rate is considerably lower (e.g., Zeebe and Wolf-Gladrow, 2001), which causes deviations from the expected carbonate parameters, as observed in the data for waters collected in December (Figure 6). The large extent of scattering for the low-salinity waters in the pH-DO% and ar-DO% plots may also reflect differences in gas exchange rate between the oxygen and CO 2 . Variability in the TA of freshwater end-members is another possible reason for the scattering.
Furthermore, the dissolution of CaCO 3 has been suggested by the high TA relative to the S values in low-bottom waters of some coastal regions (e.g., Cross et al., 2013;Wang et al., 2017). In Tokyo Bay, elevated TA values were observed ( Figure 4A) in cold waters during winter, and not in the low-ar bottom waters during the summer. This elevation from the regression line in Figure 4A was observed every year, by an average of 24 µmol kg −1 (n = 65) and up to 81 µmol kg −1 for waters at T < 13 • C. Since the ar values of these waters are high (2.8 ± 0.5) and that sedimentary dissolution of CaCO 3 is a slow process (time required for complete dissolution of available aragonite solid is 20,000 days at = 0.5, Cross et al., 2013), the elevated TA values are probably not caused by CaCO 3 dissolution. An increase in the TA during winter is not evident in offshore subtropical source water (Ishii et al., 2011) and in river water (Bureau of Waterworks Tokyo Metropolitan Government). Since surface water at higher latitudes tends to have higher TA at the same salinity level (Takatani et al., 2014), the observed increase in TA could be explained if water of northern-origin flows into Tokyo Bay every winter. However, to the best of the authors' knowledge, this has not been reported. Particulate and dissolved organic matter (POM and DOM) are known to affect the TA of coastal water samples (Cai et al., 1998;Kim et al., 2006;Hernández-Ayon et al., 2007;Kim and Lee, 2009). However, the POM and DOM concentrations observed at our study sites are high in summer and low in winter . Therefore, the elevated TA during winter cannot be explained by the organic contributions. Anaerobic remineralization processes such as denitrification and sulfate reduction can also increase the TA Krumins et al., 2013). Although increase in TA was not evident in the hypoxic/anoxic summer bottom water, sediment denitrification can be stimulated in winter by supplying O 2 to promote nitrification or NO 3 itself into nitratelimited sediment (Bonaglia et al., 2014). Quantitative analysis of the TA efflux from the sediment may clarify the reason for the TA elevation in winter waters and its importance in the carbonate system of Tokyo Bay.

Future Progress of Ocean Acidification
In this section, the future changes in the pH and ar values of surface and bottom waters at station F3 are estimated. As per the report by Yamamoto-Kawai et al. (2015), the future DIC values are calculated as DIC = DIC EQ + DIC, where DIC EQ is the DIC of the water that is in equilibrium with atmospheric CO 2 , and DIC represents an air-sea disequilibrium caused by an insufficient air-sea gas exchange, biological activities, etc. Using the CO2SYS program, the mean DIC EQ of the study period with atmospheric xCO 2 of 400 ppm is estimated from the monthly mean T obs , S obs , and TA phys . Then, DIC is calculated by subtracting DIC EQ from the monthly mean of the observed DIC. This is similar to the calculations for the carbonate parameters in the physical processes-only case and the values, with the only difference being that monthly mean values were used instead of each individual observation. Under the assumption that the DIC remains unchanged over time, DIC is added to DIC EQ estimated for the water under future xCO 2 conditions. Future changes in the pH and ar were then estimated using the derived DIC along with T obs , S obs , and TA phys values. Figure 7 shows changes in pH and ar with increasing atmospheric xCO 2 . It is predicted that in future, the decrease in the pH of surface waters will be larger in winter than that in summer ( Figure 7C). In contrast, the decrease in the pH of bottom waters will be larger in summer than that in winter, due to the lower buffer capacity of waters with higher pCO 2 . Aragonite undersaturation in bottom waters will be a common event in August when the atmospheric xCO 2 values reach 500 ppm (in the year 2043 and 2057 in RCP 8.5 and RCP6.0 scenario, respectively). When the atmospheric xCO 2 is 700 ppm (in the year 2073 in RCP8.5, after the end of this century in RCP6.0), aragonite undersaturation in the bottom waters lasts for ∼3 months from July to September (Figure 7B). Yamamoto-Kawai et al. (2015) have estimated that since 1950/1960s, eutrophication has increased the summer bottom water DIC by 46 µmol kg −1 and decreased summer surface water DIC by 46 µmol kg −1 . When we remove these changes in DIC in the calculation of future pH and ar for summer months (from June to September), results change significantly (Figure 7). Under the atmospheric xCO 2 of 700 ppm, removal of eutrophication effect decreases pH by ∼0.10 and ar by ∼0.48 in the surface water (red filled circles and red open circles). In the bottom water, the pH and ar values increase by ∼0.14 and 0. 30, respectively (red filled squares and red open squares) And aragonite undersaturation occurs only in August. In addition, in the absence of eutrophication, the rate of decrease in pH (for every 100 ppm increase in xCO 2 ) is reduced from 0.076 to 0.067 in December and from 0.102 to 0.096 in August (mean rate for xCO 2 range from 400 to 700 ppm; Figure 7). The mean ar in August becomes < 1 when the xCO 2 value reaches 650 ppm (in the year 2067 and 2095 in the RCP 8.5 and 6.0 scenarios, respectively), indicating that de-eutrophication can postpone the timing of bottom aragonite undersaturation by 24 and 38 years, RCP 8.5 and RCP6.0 scenarios, respectively. Nevertheless, deeutrophication should be accompanied with an effort to maintain low atmospheric xCO 2 values (below 650 ppm) to ensure that the water in Tokyo Bay is oversaturated with respect to aragonite throughout the year.

SUMMARY AND CONCLUSION
Seven years of observations have revealed the seasonal and interannual variations in the carbonate parameters and the processes affecting the carbonate system in Tokyo Bay. In the surface layer, the ar increases by ∼1.0 from winter to summer because of warming. The sporadic extreme freshening caused by river water input lowers the ar by ∼1.0. These physical processes do not cause large variations in pCO 2 and pH when a CO 2 equilibrium is achieved by an air-sea gas exchange. In contrast, the active photosynthetic processes in this highly eutrophic bay significantly impact the carbonate parameters. Surface water was undersaturated with respect to atmospheric CO 2 for most months, making the bay a strong CO 2 sink . Quantitative evaluation shows that photosynthesis, together with other non-physical processes, decreases the pCO 2 in surface waters by up to 345 µatm, causing both the pH and ar of surface waters to increase by up to 0.6 and 3.3, respectively. In the bottom layer, remineralization is the major process controlling the carbonate parameters. Aragonite undersaturation was observed twice during the study period in anoxic waters. The relationship between pCO 2 and DO% suggests that anaerobic remineralization under anoxic conditions lowered the ar value, causing undersaturation. This suggests that anoxic conditions can enhance ocean acidification at bottom water which can be critical to the benthic marine organisms. In addition, advection of the anoxic and aragonite undersaturated seawater to the upper water column can be detrimental to the planktonic organisms in Tokyo Bay. If the concentration of atmospheric CO 2 continues to increase, aragonite undersaturation of the summer bottom water will be a common event in the Tokyo Bay in the near future. De-eutrophication, resulting from pollution control, can decelerate ocean acidification and postpone the occurrence of prolonged aragonite undersaturation in the bottom waters. De-eutrophication will also improve the oxygen conditions in summer bottom waters, thereby suppressing anoxic remineralization and sporadic aragonite undersaturation. However, if atmospheric CO 2 exceeds 650 ppm, aragonite undersaturation will be a common feature in the summer bottom waters, even if hypoxia/anoxia are alleviated by deeutrophication.

DATA AVAILABILITY STATEMENT
The raw data supporting the conclusions of this article will be made available by the authors, without undue reservation.

AUTHOR CONTRIBUTIONS
MY-K and HK contributed the concept and design of the study. JK, MY-K, and SI collected the samples during field samplings and analyzed the samples. MY-K drafted the manuscript. All authors contributed to the manuscript revision, read, and approved the submitted version.