{'fraction_non_alphanumeric': 0.029992539169360854, 'fraction_numerical': 0.020790848047749316, 'mean_word_length': 4.702212138400454, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Dissolved oxygen variability of Indonesian seas over decades as detected by satellite remote sensing

The dissolved oxygen (DO) decrease in the ocean is a notable issue because of its potential impacts on marine biogeochemical cycles and ecosystem services. Satellite remote sensing application to support in-situ measurement is a time and cost-saving on wide scales DO monitoring. This study aims to determine the DO variability from 1993 to 2020, identify the potential areas to experience deoxygenation, and investigate the correlation between DO and other ocean parameters in Indonesian seas. The validation between in-situ and satellite-derived DO shows the determination coefficient of 0.73, indicating the satellite dataset reliability for the entire analysis. The multiple regression analysis among the long-term satellite-derived ocean parameters shows that the in-situ DO can be estimated by the combination of the potential temperature, total chlorophyll-a, and salinity. The potential temperature was statistically identified as the parameter with the highest correlation and influence on DO. The results of DO variability analysis show the overall decreasing trend with significant decreases in 1998, 2010, and 2016. There is a distinct difference in DO’s seasonal patterns in the southwestern and northeastern regions. The potential of ocean deoxygenation is detected in western Sumatra waters and the Arafura Sea at the 200–1,000 meters depth.


Introduction
Global observations conclude that the ocean lost up to 2% of dissolved oxygen (DO) over the past century [1,2,3]. Ocean deoxygenation becomes one of the most significant changes occurring in the ocean because of its impacts on biogeochemical cycles and marine ecosystem services provided by the ocean and coastal waters [4,5,6]. Global warming is the main contributor to ocean deoxygenation phenomena because the solubility of oxygen has been proven to decrease with the warming of waters [7]. Climate change also indirectly affects surface water oxygen deficiencies by reducing the ocean ventilation and extending the stratification period [8,9].
Studies of DO variability, factors influencing to the alteration of biogeochemical cycles, and the interpretation of the implications between variabilities have been crucial across disciplines these days to predict the magnitude and spatial patterns of ocean deoxygenation. However, studies at large spatial and temporal scales are needed to support the sustainability of ecosystem services provided by the ocean [10]. Current models underestimate the long-term variability of ocean deoxygenation and generally simulate only recent ocean oxygen loss inferred from few in-situ observations [5]. Recently, the application of satellite earth observation data is more frequently applied in water quality studies [11,12,13]. The advantages of remotely sensed analysis over in-situ measurement are the extended temporal frequency and the greater spatial coverage with better resolution [7,13]. This method will be beneficial, effective, and efficient for the broader study area.
Indonesian seas represent the climate-regulated marine region, especially with the influence of upwelling systems and ocean water exchanges in several areas. This region will be a suitable and engaging locus study to understand the variabilities of dissolved oxygen. This study aims to understand DO variability from 1993 to 2020, identify the potential areas to experience deoxygenation, and investigate the correlation between DO and other ocean parameters in Indonesian seas using in-situ and satellite data. The outputs resulted from this study are expected to provide a database and scientific information to be utilized in the development of climate change mitigation and adaptation measures. This study will also show the potential capability of the remote sensing application in ocean deoxygenation monitoring both in high spatial and temporal resolution.

Study Site
The Indonesian territorial waters are located in tropical latitudes between two continents (Asia and Australia) and two oceans (Pacific and Indian). As the locus study, Indonesian seas are divided into 4 (four) regions ( Figure 1). Region A is along the coast of Sumatra, Java, Bali, and Lombok islands that directly connected to the Indian Ocean, located in 91.3⁰ E -107.2⁰ E and 7.3⁰ N -11.0⁰ S. Region B is an area in the southern part of eastern Indonesian waters, located in 117.2⁰ E -141.0⁰ E and 1.5⁰ S -8.0⁰ S. Region C is an area in the northern part of eastern Indonesian waters directly connecting to the Pacific Ocean, located in 117.2⁰ E -141.0⁰ E and 5.0⁰ N -1.5⁰ S. Lastly, Region D is an area in the Sunda Shelf including in Indonesian territorial, located in 95.0⁰ E -107.2⁰ E and 7.3⁰ N -8.5⁰ S.

Satellite Data
Four variables used in this study were derived from satellite data, i.e., dissolved oxygen (DO), total chlorophyll-a (CHL-a), potential temperature (PT), and salinity (SAL). DO and CHL-a data are retrieved from Global Ocean Biogeochemistry Hindcast and Global Ocean Biogeochemistry Analysis and Forecast. PT and SAL data are retrieved from Global Ocean Physics Reanalysis and Global Ocean 1/12° Physics Analysis and Forecast Updated Daily. These Level 4 data are provided by E.U. Copernicus Marine Service Information from 1993 to 2020 with the spatial resolutions of 0.25 degrees for biogeochemistry data and 0.083 degrees for physics data and vertical coverage from 0 -5,500 m below sea level, except for Region D that only covered up to 50 m below sea level.  Figure 1. Region A, B, C, and D of Indonesian seas as the study area.

In-situ Data
In-situ data derived from World Ocean Database ver.18 (WOD18) provided by National Centers for Environmental Information (NCEI), National Oceanic and Atmospheric Administration (NOAA). Unfortunately, no in-situ data was available in Region D. The variables and periods of WOD18 data used are in similar kinds and ranges with the hindcast satellite data. DO, PT and SAL data were retrieved by Conductivity-Temperature-Depth (CTD) instrument and CHL-a data retrieved by Ocean Station Data (OSD) instrument. Each data value and each profile in WOD18 has quality control flags associated with the automatic checks system.

Data Analysis
To evaluate quantitatively the spatial-temporal variability of DO and produce the long-term trends of the ocean deoxygenation in Indonesian seas, the monthly average values of satellite data in specified depth were extracted and visualized in charts and maps. The vulnerable deoxygenation areas, Oxygen Minimum Zone (OMZ), and extreme events of ocean deoxygenation were identified from the charts and the maps. All satellite data processing were processed by QGIS ver.3.16 software, while the in-situ data were processed by Ocean Data View ver.5.4.0 software. Statistical analyses were performed by R ver.3.6.1 software.
The in-situ data of DO were compared with remote-sensing estimates in the same periods, locations, and depths. Since, the in-situ data and shallow covers of satellite data in Region D was unavailable, the analysis was only based on the sea surface. Linear regression analysis between the matching values of satellite and in-situ DO data were used to assess the validation of satellite data and its application to the spatial-temporal analyses. A matching value was based on the closest 3 x 3 satellite pixel values to the in-situ measurement location within a one-day time interval. The 6 of 9 satellite pixels were used as a minimal validation and the average difference between the central pixel and all other pixels was less than 25%. The coefficient determination (R 2 ) was used as significant relationship in the linear regression statistical analysis to use the reliability of the entire analysis. Mean Absolute Deviation (MAD), Root Mean Square Error (RMSE), and Mean Absolute Percentage Error (MAPE) were also calculated to find out the error degree of satellite data in predicting in-situ value. After the surface and subsurface analyses were executed, the Pearson Product-Moment Correlation analysis was used to find the correlation value between DO and the other three variables. The Multiple Linear Regression analysis was also performed to find the most influential variable on DO variability and to develop estimated DO algorithm in each region.

The reliability of the satellite data
The matching value procedures produced 220 matching pairs of satellite and in-situ DO concentrations of Indonesian seas. The validation of these matching pairs shows the coefficient determination (R 2 ) of 0.73, MAD of 1.92, RMSE of 2.10, and MAPE of 0.95%. Regarding to the remote sensing-based models experimental for water quality parameters [14], these satellite data are quite reliable to be used in the entire analysis. Figure 2 shows the comparison of in-situ and satellite-derived DO concentration in Region A, Region B, and Region C. Considering the high validation value between in-situ and satellite data of Regions A, B, and C, we assume that satellite data in Region D is also valid and feasible to be analyzed in this study.  Figure 3 shows the monthly concentration of surface DO in Region A, B, C, and D in 1993 -2020. In general, there is evidence of ocean deoxygenation in the four regions with the most distinct DO declining trend in Region A. The highest annual average of surface DO concentration from 1993 -2020 found in Region A, followed slightly by Region D, then Region C and Region B, by 202.34 ± 2. 22 Significant seasonal patterns are shown in Region A and B that are connected directly to the Indian Ocean monsoonal winds. DO concentration in both regions increases during the west monsoon (wet season), with the highest peak occurs in October. During the east monsoon (dry season), DO concentration will decrease with the lowest point in April. Meanwhile, Region C is showing no distinct seasonal pattern of DO and tends to be flat. The increasing and decreasing of DO occur twice a year (high in April and October and low in June and December) with concentrations that tend to be flat. On the other hand, Region D shows another different seasonal pattern, which is the highest DO concentration found in February (wet season) and lowest in August (dry season), similar to the seasonal influence in Region A and B but with different monsoonal pattern. The DO seasonal pattern of Region D is similar to those in higher latitude regions such as the Arabian Sea [17] and Yellow Sea [7].  Figure 3). According to the El Niño Southern Oscillation Index (SOI) and Indian Ocean Dipole (IOD) index [15,16], these decreases occurred during strong ENSO strengthened by association with both positive and negative IOD. Strong ENSO values were recorded in Indonesian waters from 1993 to 2020, with El Niño and positive IOD events occurred in 1997 -1998 and 2015 -2016, and La Niña and negative IOD events occurred in 2010 -2011. The causative correlation of DO decreases and ENSO-IOD events will not be explained in this study since there is no available direct data to support these co-occurrences. However, based on the coincidences between the three extreme events recorded with the significant DO decreases mentioned previously, it can be assumed that the extreme climate phenomena will affect DO concentration, either directly or indirectly, on the sea surface.

Extension of oxygen minimum zone indicates the ocean deoxygenation in the Indonesian seas
The lateral and vertical distribution of DO can be detected from the depth section of DO level ( Figure  4). Based on these depth sections, we identified the temporal changes in DO concentration and recognized the locations of the most significant ocean deoxygenation. The most significant deoxygenation pattern was found in Region B, which centered around the Arafura Sea in southern Papua Island. In 1993, the lowest DO concentration in the water columns of the Arafura Sea was above 86 mmol/m 3 . However, the DO concentration found dropped drastically to 0 -20 mmol/m 3 in 2020, at the depths of 240 -1,100 m below sea level and extends westward. The water column transforms to OMZ with a deoxygenation rate of 2.6 mmol/m 3 /year. In comparison, this low concentration of DO is similar to the concentrations in the Arabian Sea and Bay of Bengal, with 0 -20 mmol/m 3 located at the middepths of the water column [17,18].  The deoxygenation area in Region A is located in the western Sumatra waters with a deoxygenation pattern that extends to the southern Java Island. In the northernmost part of Sumatra waters, the OMZ extends laterally eastward with a deoxygenation rate of 1.3 mmol/m 3 /year and reached its lowest point of 0 -20 mmol/m 3 in 2020 at the depths of 180 -700 m below sea level. Compared to Region A and Region B, Region C seems not experienced deoxygenation. The DO concentration actually increased in 2020 to 140 mmol/m 3 in the depth of 600 m below sea level, meanwhile in 1993, the concentration on the depth only reached 110 mmol/m 3 .
We collected CTD data of DO, PT, and SAL of Indonesian seas from the surface to 3,000 m depth below sea level. Figure 5 shows the comparison of DO with PT (left) and SAL (right) along with the depths. The highest DO concentration is always found between 0 -150 m depth. The two main reasons why oxygen content is highest at the sea surface are because it is the location where oxygen dissolves into the ocean from the atmosphere, and it is also the location where oxygen is produced by phytoplankton through photosynthesis [19]. Detected deoxygenation in Indonesian seas dominates the depth of 200 -1,300 m below sea level, with situations of higher temperature and salinity.

Significant correlation of dissolve oxygen, temperature, and salinity
The oxygen deficiency of the Indonesian sea surface most likely occurs in the east monsoon (dry season) associated with the increases in water temperature and salinity. On the contrary, the temperature will decline in the west monsoon (wet season) and affect the enhances of dissolved oxygen on the seas. Another notable factor is the riverine inputs that are higher in the rainy season. In the west monsoon, the heavy and intense rains will be cooling down the temperature, and the flows will bring out the saline level towards the open sea. The terrestrial-sourced chlorophyll-a and nutrient concentrations are also carried away by the riverine discharge to the estuary and open ocean. This phenomenon explains why the seasonal pattern of the surface DO in the Indonesian seas has an opposite pattern with the surface PT and SAL, whilst the similar pattern is noticeable in surface CHL ( Figure 6). As the PT and SAL decrease, it will affect the ability of organisms to obtain oxygen from the water. PT and SAL independently affect the solubility of seawater, which is also controlled by the partial pressure of the oxygen (PO2) [20]. The increase in PO2 will enhance solubility and diffusion of the seawater and will slow the rate of oxygen deficiency [21].   differences. In these two regions, PT and DO have dominant strong negative correlations (Pearson correlations close to -1.00), indicating inverse relations. These relationships mean that if the sea temperature of the regions increases, the oxygen concentration in the ocean will decrease, supporting the theory of ocean deoxygenation worsened by global warming. However, an anomaly was found at the depths of 50 m below sea level, when the sea temperature gives positive correlations to the dissolved oxygen. As stated by [22], this kind of anomaly may occur in the tropical seawater at the depth 50 -75 m, due to the dynamic process of heating-cooling of temperature combined with production-entrainment of oxygen. The other main driver also wind-driven adiabatic displacements of the thermocline.
Such inverse relationships are also noticeable between DO and SAL in Region A and Region B and stronger with the depth. These negative correlations mean that increasing the salt concentration in the ocean leads to oxygen content decreases, and salinity level has a more crucial impact on oxygen solubility in the deep sea than in the shallow sea. Chlorophyll-a was more concentrated on the sea surface where higher oxygen concentration and higher nutrient levels were supplied. The correlation between DO and CHL-a in the two regions gives not constant values and becomes weaker as depth increases. Strong positive correlations are found in the sea surface, become negative in the depths of 50 -100 m, and change into positive again with weak impacts (Pearson correlations between 0 -0.50).
Region C shows a different correlation pattern with the previous two regions. DO has a negative correlation with PT at the depths of 0 -50 m below sea level but strong only at the sea surface (Pearson correlation -0.60). From the depth of 100 -1,000 m below sea level, DO and PT relationship change to a weak positive correlation. Meanwhile, salinity level always gives negative correlations to oxygen solubility at all depths but are weak. In contrast to Region A and Region B, the chlorophyll-a concentration at sea level Region C almost provides no correlation to oxygen solubility (Pearson correlation -0.06). Then, the correlations become positive at the depths of 100 -200 m below sea level, and changes become negative correlation to a depth of 1,000 m below sea level. We assumed that the variation in DO, PT, SAL, and CHL-a correlations in these regions was affected by external factors, including weather, wind, and currents.  How much the PT, SAL, and CHL-a variables were able to explain the DO variation were analyzed by calculating the coefficient of determination (R 2 ). By tested the in-situ data found in the Indonesian seas, the results show that PT individually has the largest influence on DO at 66.4%, followed by SAL and CHL-a at 46.8% and 30.4%, respectively. However, the effect of these three variables on DO becomes greater when they are combined. The combination of PT and SAL affects DO by 76.0%, much greater than the individual effect of the two variables. Moreover, the combination of PT, SAL, and CHLa affects DO by 76.3%, which means that there are other effects of 23.7% from other variables that not explained in this study. The test proves that PT, SAL, and CHL-a variables contribute significantly to explaining DO variations in the Indonesian seas. As reported from the previous study [23,24], that PT, SAL, and CHL-a will eventually change over time, we can expect the significance change of DO in the future that may cause the vast impact of deoxygenation.

Multiple linear regression analysis of DO with satellite-derived variables
Three satellite-derived variables of PT, SAL, and CHL-a were used to estimate DO in Region A, B, C, and D by applying multiple linear regression analyses. By inputting the variables with the stepwise regression method, the models with the highest determination coefficient and free of multicollinearity issues were developed. However, we could only estimate the surface DO in Region D due to the limited vertical data in this area. The results show strong correlations of predicted in-situ DO with linear combinations of satellite-derived PT, SAL, and CHL-a for Region A and B, and linear combinations of satellite-derived PT and SAL for Region C and D (Table 1). Variable CHL-a in the model of Region C and D was removed because its p-value is far larger than 0.05 and statistically assumed to have no significant contribution on DO (p-value = 0.65 and 0.17 in Region C and D, respectively). Based on the results of the multiple linear regression analysis, an empirical DO algorithm can be derived as: Where DOA is estimated DO in Region A, DOB is estimated DO in Region B, DOC is estimated DO in Region C, and DODS is estimated surface DO in Region D. Using the satellite-derived DO algorithm, long-term changes in DO concentration can be detected in the areas with limited observation data. Nevertheless, in-situ data from Region D is still needed for validation purpose in future studies.