Automatic Mapping and Monitoring of Marine Water Quality Parameters in Hong Kong Using Sentinel-2 Image Time-Series and Google Earth Engine Cloud Computing

Kwong, Ivan H. Y.; Wong, Frankie K. K.; Fung, Tung

doi:10.3389/fmars.2022.871470

ORIGINAL RESEARCH article

Front. Mar. Sci., 10 May 2022

Sec. Ocean Observation

Volume 9 - 2022 | https://doi.org/10.3389/fmars.2022.871470

Automatic Mapping and Monitoring of Marine Water Quality Parameters in Hong Kong Using Sentinel-2 Image Time-Series and Google Earth Engine Cloud Computing

Ivan H. Y. Kwong^1,2*

Frankie K. K. Wong¹

Tung Fung^1,2

¹Department of Geography and Resource Management, The Chinese University of Hong Kong, Hong Kong, Hong Kong SAR, China
²Institute of Future Cities, The Chinese University of Hong Kong, Hong Kong, Hong Kong SAR, China

Continuous monitoring of coastal water qualities is critical for water resource management and marine ecosystem sustainability. While remote sensing data such as Sentinel-2 satellite imagery routinely provide high-resolution observations for time-series analysis, the cloud-based Google Earth Engine (GEE) platform supports simple image retrieval and large-scale processing. Using coastal waters of Hong Kong as the study area, this study utilized GEE to (i) query and pre-process all Sentinel-2 observations that coincided with in situ measurements; (ii) extract the spectra to develop empirical models for water quality parameters using artificial neural networks; and (iii) visualize the results using spatial distribution maps, time-series charts and an online application. The modeling workflow was applied to 22 water quality parameters and the results suggested the potential to predict the levels of several nutrients and inorganic constituents. In-depth analyses were conducted for chlorophyll-a, suspended solids and turbidity which produced high correlations between the predicted and observed values when validated with an independent dataset. The selected input variables followed spectral characteristics of the optical constituents. The results were considered more robust compared to previous works in the same region due to the automatic extraction of all available images and larger number of observations from different years and months. Besides visualizing long-term spatial and temporal variabilities through distribution maps and time-series charts, potential anomalies in the monitoring period including algal bloom could also be captured using the models developed from historical data. An online application was created to allow novice users to explore and analyze water quality trends with a simple web interface. The integrated use of remotely-sensed images, in situ measurements and cloud computing can offer new opportunities for implementing effective monitoring programs and understanding water quality dynamics. Although the obtained levels of accuracies were below the desired standard, the end-to-end cloud computing workflow demonstrated in this study should be further investigated considering the cost and computational efficiency for timely information delivery.

Introduction

Water quality has gradually experienced degradation in many coastal regions because of pollution and heavy use resulting from anthropogenic activities (Chen et al., 2004; Gholizadeh et al., 2016). In Hong Kong, the Harbour Area Treatment Scheme has been commenced since 2001 to improve the water quality of Victoria Harbour at a cost of US$ 2.5 billion (Xu et al., 2011). Nevertheless, the marine environment of Hong Kong is still under considerable pressure from local, regional and global stressors such as reclamation, overfishing, biological invasion, transboundary pollution and climate change (Lai et al., 2016). Continuous monitoring and assessment of the water quality are critical for water resource management and marine ecosystem sustainability.

Water quality can be considered as a general descriptor of water properties in terms of physical, chemical and biological characteristics (Ritchie et al., 2003). Many water quality parameters, such as chlorophyll-a (Chl-a), suspended solids (SS), and turbidity, have been commonly used as indicators (Hafeez et al., 2019; Topp et al., 2020). For instance, Chl-a is a photosynthetic pigment and can be used as indirect measurements of phytoplankton biomass, algal production and trophic state (Schalles, 2006; Flores-Anderson et al., 2020). Rapid growth of algae or cyanobacteria can lead to harmful algal blooms that adversely affect people and ecology (Poddar et al., 2019; Wang et al., 2020). SS refers to both inorganic and organic particles held in suspension throughout a water column (Topp et al., 2020). Their movements are related to sediment transport, nutrient cycle and water clarity and are often linked to economic activities in ports and waterways (Shahzad et al., 2018; Wang et al., 2018). Turbidity is closely related to SS and is a measurement of light scattering within a water column caused by particles (Matthews, 2011). Increasing turbidity can reduce light penetration and cause cloudiness of the water (Khan et al., 2021).

Traditionally, water samples are collected from the field and then analyzed in the laboratory. Although such in situ measurements offer high accuracy, they are costly, labor-intensive and time-consuming (Gholizadeh et al., 2016; Pizani et al., 2020). Spatial and temporal variations within the entire waterbody are also hard to be identified from the sampling points (Ouma et al., 2020). On the other hand, recent advances in remote sensing technology have enabled the extraction of spatial and temporal profiles of surface water by providing a synoptic and synchronized view from aerial or satellite platforms over a large area (Schaeffer et al., 2013; Ansper and Alikas, 2019; Flores-Anderson et al., 2020). Since the upwelling radiance from the water at various wavelengths often consists of signals from the optically active components of the water (Ritchie et al., 2003; Shahzad et al., 2018), the spectral shapes and magnitudes of the water-leaving reflectance can be used directly or indirectly to detect different water quality parameters through empirical or analytical models (Matthews, 2011; Salem et al., 2017; Wang and Yang, 2019).

Despite a large amount of satellite data being available for remote sensing of water quality, only a few satellites can provide sufficient spatial, spectral and temporal resolutions to resolve the complex water variability in coastal areas. Ocean color sensors specifically designed for marine purposes (e.g., SeaWiFS, MODIS, MERIS) (Hu et al., 2012) cannot fulfil the requirements of coastal water monitoring due to their coarse spatial resolutions (> 250 m) (Mouw et al., 2015). Land sensors with moderate resolutions, especially the Landsat series satellites (30 m), were found to be the most popular data resources, mainly because of their temporal coverage and easy accessibility (Wang and Yang, 2019). Still, the revisit time of Landsat-8 (16 days) limits its use in routine monitoring of water quality (Toming et al., 2016).

Currently, the Multispectral Instrument (MSI) onboard Sentinel-2, launched by the European Space Agency Copernicus program in 2015, has opened a new potential in water remote sensing (Drusch et al., 2012). The MSI sensor measures radiance in 13 spectral bands spanning from visible to shortwave infrared with 10–60 m spatial resolution, including three red-edge bands that can be advantageous over optically complex coastal waters (Pahlevan et al., 2017; Ansper and Alikas, 2019), and provides a 5-day revisit frequency with both Sentinel-2A and Sentinel-2B satellites (Toming et al., 2016). In particular geographical areas like Hong Kong which are located at the edge of the sensor swath, Sentinel-2 is also preferred over Landsat satellites since its wider swath width can provide complete coverage of the study area in a single scene.

Images obtained from Sentinel-2 MSI have been demonstrated to successfully derive water quality parameters in a variety of scenarios, such as inland water reservoirs (Ouma et al., 2020; Pizani et al., 2020; Pompêo et al., 2021), lakes (Liu et al., 2017; Ansper and Alikas, 2019; Free et al., 2020; Soomets et al., 2020), rivers (Kuhn et al., 2019), estuaries (Sent et al., 2021) and bays with aquaculture (Gernez et al., 2017; Soriano-González et al., 2019). However, most of these analyses focused on the development of algorithms and only utilized images acquired from a single date or a limited temporal range (Topp et al., 2020). The strength of remote sensing lies in providing both comprehensive historical records of water quality, which facilitates understanding of long-term environmental changes (Kim et al., 2020), and near real-time measurement at critical points, where in situ measurements are not readily available (Sobel et al., 2020). Considering that Sentinel-2 data will be available routinely for many years and free of charge, several studies also emphasized the importance of using multi-temporal images to provide time-series analysis of changes and trends (Toming et al., 2016; Tamiminia et al., 2020; Elhag et al., 2021).

When utilizing satellite images for routine monitoring of water quality parameters, researchers have argued that computational demands required for managing and analyzing large volumes of data can be a major bottleneck (Kumar and Mutanga, 2018; Page et al., 2019). Compared to conventional image processing techniques, Google Earth Engine (GEE), a cloud-based platform for large-scale geospatial analysis (Gorelick et al., 2017), has been demonstrated to provide easy accessibility to readily available data and powerful computing capabilities for a variety of applications (Traganos et al., 2018; Tamiminia et al., 2020; Jia et al., 2021; Li et al., 2021). However, only a few studies have attempted to exploit the functionality of GEE to develop an automated approach or even a web application that can further enhance usability of the outputs (Markert et al., 2018; Rudiyanto et al., 2019; Canty et al., 2020). Hence, this study aimed to fill the gap by developing a low-cost, end-to-end and efficient workflow that can provide rapid information for improving water quality monitoring.

Among the previous remote sensing-based studies of water quality in Hong Kong, Wong et al. (2007; 2008) used MODIS images to retrieve Chl-a, SS and turbidity for coastal waters of Hong Kong, but the coarse spatial resolution (500 m) limited their usability in the spatially complex environment. Tian et al. (2014) produced the first map of SS concentration at moderate spatial resolution (30 m) in the Deep Bay area of Hong Kong using HJ-1 satellite images, despite the fact that their study focused only on a narrow range of SS concentrations and a specific area inside Hong Kong. Nazeer and Nichol (2015) developed an empirical model for estimation of SS concentrations by combining Landsat-5, Landsat-7 and HJ-1 images, although they are old-generation satellites with less desirable band designations and data bit levels. A recent study by Hafeez et al. (2019) focused on the comparison of several machine learning techniques for improving water quality estimation with Landsat series satellites, while Hafeez and Wong (2019) tested the performance of using Sentinel-2 images to estimate Chl-a and SS in a single year. These studies were either limited to using few image dates and in situ samples for model development and validation, or not providing a solution to operational water monitoring in different scenarios.

The objective of this study was to develop an automated framework to retrieve the spatial and temporal patterns of different water quality parameters, specifically Chl-a, SS and turbidity, for coastal waters of Hong Kong by integrating Sentinel-2 image time-series, in situ measurements from monitoring stations and the GEE cloud computing capability. In particular, this study aimed to utilize GEE to (i) query and pre-process all Sentinel-2 observations that coincided with in situ measurements; (ii) extract the spectra to develop empirical models for water quality parameters using artificial neural networks (ANN); and (iii) visualize the results using spatial distribution maps, time-series charts and an online application. The first part of this study presents the processing workflow and demonstrates that it could be applied to all 22 water quality parameters collected by monitoring stations. The second part of this study focuses on the modeling results of three major water quality parameters, including Chl-a, SS and turbidity. Through evaluating the empirical modeling performance and analyzing the water quality changes, this study provided both theoretical and practical implications, facilitating the implementation of an effective water resource monitoring program.

Materials and Methods

Study Area

Hong Kong is a highly urbanized city lying at the northern limits of the Asian tropics between latitudes 22°08’N and 22°35’N and longitudes 113°49’E and 114°31’E (Figure 1). The climate is subtropical, with hot rainy summers from May to August and cool dry winters from November to March. Situated with the Pearl River Estuary to the west and surrounded by the South China Sea to the south and east, Hong Kong has a land area of 1,108 km² including more than 250 offshore islands. The coastline extends around 1,200 km and the sea area is 1,648 km². The marine waters support diverse forms of marine life ranging from microscopic algae to dolphins, and are used for navigation, recreation, seafood production, and supply of flushing and cooling water (Environmental Protection Department, 2021).

FIGURE 1

Figure 1 (A) The study area of this study, including the 10 water control zones defined by the Environmental Protection Department (EPD) of the Hong Kong government; (B) Location of Hong Kong in the regional context.

Due to the combination of anthropogenic activities and variable hydrographical conditions, the coastal waters of Hong Kong are physically and chemically complex (Hafeez et al., 2019). For example, in the western areas (e.g. Deep Bay and North Western Zones), estuarine waters are influenced by freshwater discharge from the Pearl River. Water quality in this area is also influenced by the sediment loads and nutrients in the river discharges (Zhou et al., 2007). Conversely, the eastern waters (e.g. Mirs Bay) are more oceanic and influenced by the Pacific currents (Lai et al., 2016). The central areas represent a transition zone that varies seasonally, in addition to the influence of local pollution load from Victoria Harbour (Xu et al., 2011). The water quality also varies spatiotemporally depending on two seasonal monsoons. While the northeast monsoon prevails and enhances the effect of China Coastal Current during winter, the southwest monsoon in summer transports continental shelf water landward (Zhou et al., 2012). The estuarine plume covers the northern part of South China Sea when the river discharge is large in summer. The normal tidal range in Hong Kong waters is between one and two meters, with flood tide for oceanic water to flow to the northwestern direction and ebb tide in the reversed direction.

Water Quality Parameters and Station-Based Data

The Environmental Protection Department (EPD) of the Hong Kong government has divided Hong Kong waters into 10 water control zones (Figure 1) based on the hydrodynamic characteristics and pollution status. A systematic marine water quality monitoring program has been in operation since 1986 to measure a total of 22 water quality parameters (Environmental Protection Department, 2021). The list of the measured parameters is shown in Table 1, and they are related to different aspects of water quality including (i) physical and aggregate properties, (ii) aggregate organic constituents, (iii) nutrients and inorganic constituents, and (iv) biological and microbiological examination. The modeling workflow was applied to all 22 parameters in the first part of the present study, while detailed analyses focused on three major parameters, including Chl-a, SS and turbidity.

TABLE 1

Table 1 Summary statistics of station-based monitoring data of the water quality parameters extracted in this study.

EPD takes measurements and collects water samples at 76 monitoring stations in open waters every month using their dedicated marine monitoring vessel equipped with Differential Global Positioning Systems (DGPS) and an advanced conductivity-temperature-depth (CTD) profiler. The water quality parameters are measured from three depths, including surface water (1 m below the sea surface), middle water (half the sea depth), and bottom water (1 m above the seabed). Water samples are collected in a 500 mL Nalgene bottle and analyzed by EPD’s in-house laboratory. To measure Chl-a concentration (µg/L), the American Public Health Association (APHA) 20ed 10200H2 spectrophotometric method based in-house GL-OR-34 method is used. For SS concentration (mg/L), the APHA 22ed 2540D weighing method based in-house GL-PH-23 method is used, while turbidity (NTU) is measured on-site by the OBS-3 turbidity sensor linked to a SEACAT 19+ CTD and Water Quality Profiler.

The marine water quality monitoring data are open data in the Hong Kong government online database (https://data.gov.hk/). The data are updated yearly and released usually in the third quarter of the next year. In this study, we retrieved the water quality data from 2015 to 2020 which match the availability of Sentinel-2 images. While this study was conducted in late 2021, station-based data covering the current year were not yet available. A summary of the station-based data is provided in Table 1.

Selection and Pre-Processing of Sentinel-2 Images

Sentinel-2 is an Earth observation mission from the Copernicus Program and consists of two identical satellites. Sentinel-2A was launched on 23^rd June 2015 and Sentinel-2B was launched on 7^th March 2017. The band configuration of Sentinel-2 is given in Table 2. The bands range from 443 nm to 2190 nm, featuring four bands at 10-m (Visible and near-infrared [NIR]), six bands at 20-m (red-edge and shortwave infrared [SWIR]) and three bands at 60-m (atmospheric) spatial resolution.

TABLE 2

Table 2 Spectral bands for Sentinel-2 sensors (Drusch et al., 2012). Bands 8, 9 and 10 were eliminated in this study due to the wavelength overlap or absence of water surface information.

In this study, selection of suitable Sentinel-2 images and subsequent pre-processing steps were done using GEE (https://code.earthengine.google.com/). GEE includes both a web-based Code Editor in JavaScript language and a Python application programming interface (API) for data analysis outside the web environment (Tamiminia et al., 2020). All the analysis procedures in this study were written using GEE Python API in Google Colab (https://colab.research.google.com/) such that the codes can be directly executed through a web browser with minimal configuration.

The complete archive of Sentinel-2 images (Level-1C) from 2015 is available as an image collection in the GEE repository. Despite the ready availability of atmospherically corrected Level-2A products, this collection only contains images dating back to 2017, which makes it less suitable for time-series analysis. Instead, Level-1C images covering the study area were queried in GEE with a maximum cloudy pixel percentage of 20% according to metadata. A total of 120 cloud-limited scenes from 23^rd October 2015 to 31^st December 2021, covering all seasons of the year, were retrieved and used in this study (Table 3).

TABLE 3

Table 3 Dates of the cloud-limited Sentinel-2 scenes retrieved and analyzed in this study.

The Level-1C Sentinel-2 data contained spectral values representing the top of atmosphere (TOA) reflectance. Atmospheric correction is an important process to obtain water-leaving reflectance by removing atmospheric interference from the total signal received by the satellite sensor (Soriano-González et al., 2019). This study implemented an atmospheric correction model called Py6S (Wilson, 2013), which is a Python interface of the “Second Simulation of the Satellite Signal in the Solar Spectrum” (6S) radiative transfer model (Vermote et al., 1997), to all Sentinel-2 images with the maritime aerosol setting. Py6S was chosen since it was found to be the best method over complex coastal waters of Hong Kong when tested using satellite images with similar resolutions (Nazeer et al., 2014) and relevant codes have been developed for GEE (Murphy, 2020).

For each atmospherically corrected image, a cloud mask was applied based on the s2cloudless dataset, which is a machine learning-based cloud detector precomputed on GEE (Zupanc, 2017). Cloud shadows were also estimated and masked using the position of detected clouds and the viewing geometry. Since some of the images were found to be affected by sun glints, a simple correction method was applied by subtracting 50% of the reflectance values of band 11 (SWIR) in each pixel from all bands. It was based on the assumption that the water-leaving radiance in SWIR wavelength is very low, thus a significant proportion of the SWIR signals that appear on the image could be attributed to sun glints and these signals correlate well with the glints experienced by other channels (Kay et al., 2009). The value was selected with reference to previous analysis of surface reflectance (Kuhn et al., 2019) and manual inspection of images acquired in adjacent swath with different viewing angles. In addition, the Modified Normalized Difference Water Index (MNDWI) for open waters (Xu, 2006), calculated as the difference-sum ratio of SWIR and visible green bands, was applied to all images to separate water pixels from land areas and remove low-quality pixels with a threshold value of zero.

Among the 13 spectral bands acquired by Sentinel-2, bands 9 (water vapor) and 10 (cirrus) were eliminated as they do not contain water surface information. Band 8 (NIR) was also removed since its wavelength overlaps with band 8A and the latter band provides more precise spectral measurement. The pre-processed images contained 10 spectral bands ranging from visible and NIR to SWIR wavelengths. The difference in spatial resolutions among spectral bands was not an issue in GEE which uses scale specified by the output to determine the appropriate level of input image pyramids.

Development of Empirical Models

Sentinel-2 reflectance values were extracted from the location of in situ measurements with an allowance of one-day difference, to ensure sufficient match-up points and consider the short water residence time of around two days in the wet season (Zhou et al., 2012). A 20-m buffer was adopted for each sampling point to reduce the effects of positional accuracy and noise from the sensor. All sampling locations were at least 200 m away from land to reduce adjacency effects. Furthermore, low-quality points were removed by Tukey’s fences method, which considered the point as an outlier if the value in any one of the spectral bands exceeded a distance of “1.5 times interquartile range” below the first or above the third quartiles. This resulted in a total of 352 observations, for which 300 observations in 2015–2019 were used for training and model development, while the remaining 52 observations in 2020 were used for validation, in order to demonstrate the model performance to estimate current water qualities from historical data.

Many algorithms have been proposed for retrieving water quality parameters from remote sensing reflectance. Among them, empirical methods aim to establish the relationships using statistical techniques and have advantages of computational simplicity and ease of implementation for different water quality parameters (Ritchie et al., 2003; Wang and Yang, 2019). A locally optimized algorithm is often recommended to provide better performances when calibrated with in situ observations (Tian et al., 2014; Yoon et al., 2019; Flores-Anderson et al., 2020). To establish a suitable model in this study, the following band ratios and arithmetic variables, which were developed with reference to bio-optical properties of waters (Matthews, 2011), were computed as input predictors.

Two-band ratios—The use of band ratios theoretically reduces the effects of tidal or seasonal variations and maximizes the sensitivity to water quality parameters against other constituents (Gons, 1999). For instance, a blue-green ratio uses the strong light absorption by carotenoids at 490 nm (blue) and minimal absorption of photosynthetic pigments at 560 nm (green). The ratio between 705 nm and 665 nm spectra is also based on the interaction between backscattering from particulate matter and the absorption features (Toming et al., 2016). Since there is little consistency in the optimal band ratios in different studies (Matthews, 2011), possible ratio combinations of all 10 spectral bands were computed and tested in this study. Instead of simple ratios, normalized band ratios were adopted in this study in order to limit the range between -1 and 1 (Equation 1).

\begin{array}{l} N o r m a l i z e d r a t i o (i, j) = \frac{R (i) - R (j)}{R (i) + R (j)} & (1) \end{array}

Where i and j are any bands from 1 to 12 except 9 and 10, R(i) represents the reflectance of band i.

Three-band ratio—The three-band method was developed by Dall’Olmo and Gitelson (2005) with a basic principle to find the three-band combinations that are most relevant to the absorption coefficient of the water quality parameters (Wang and Yang, 2019). A typical example is the use of 665-nm, 705-nm and 740-nm wavelengths, which showed good agreement with in situ Chl-a data in turbid and productive waters (Ansper and Alikas, 2019). In this study, all possible combinations from three consecutive bands were calculated according to Equation 2.

\begin{array}{l} T h r e e b a n d r a t i o (i, j, k) = [\frac{1}{R (i)} - \frac{1}{R (j)}] \times R (k) & (2) \end{array}

Where i, j and k are any three consecutive bands from 1 to 12 except 9 and 10, R(i) represents the reflectance of band i.

Line-height variable—The line-height algorithm is based on band differences and measures the height of a reflectance peak at a specific wavelength from a linear baseline drawn between two sides of the peak (Gitelson et al., 1994). The fluorescence line height algorithm (Gower et al., 1999) and maximum chlorophyll index (Gower et al., 2005) are examples developed using the peaks of different MERIS bands related to Chl-a concentrations. To test the applicability of various choices of peak wavelengths, line-height variables were computed in this study from each of the 10 Sentinel-2 spectral bands with the adjacent two bands (Equation 3).

\begin{array}{l} L i n e h e i g h t (i, j, k) = R (j) - R (i) - [R (k) - R (i)] \times \frac{λ (j) - λ (i)}{λ (k) - λ (i)} & (3) \end{array}

Where i, j and k are any three consecutive bands from 1 to 12 except 9 and 10, R(i) represents the reflectance of band i, λ(i) represents the central wavelength of band i.

The same empirical modeling procedures were then applied to each of the 22 water quality parameters. ANN is a supervised machine learning algorithm that consists of highly interconnected neurons to create weighted links between input data and output information (Atkinson and Tatnall, 1997). The typical model, called multilayer perceptron (MLP), uses a feed-forward connection arranged in input, hidden and output layers where summation and activation functions were performed. This type of model has been widely employed to estimate various water quality parameters (Chen et al., 2004; Chebud et al., 2012; Sagan et al., 2020) and was found to produce the highest accuracy when applied to Hong Kong waters previously (Hafeez et al., 2019). Recent studies also successfully utilized ANN models on time-series images created from GEE (Amani et al., 2020; Ghorbanian et al., 2022).

In the present study, the input layer consisted of the spectral bands and combinations, while the output layer contained the values of target water quality parameters. The candidate predictors included the original Sentinel-2 bands, their square and cubic values, as well as the three types of variables computed above. To determine an appropriate model architecture that could achieve the best possible accuracy, a 5-fold cross-validation on the training dataset was performed to select the optimal combination of parameters, including (i) five to twelve input variables which were chosen according to correlations with the dependent variable, (ii) one or two hidden layers with two to ten neurons in each layer, and (iii) L2 regularization parameter from 10^-1 to 10^-6 which penalizes large weights. The ANN models adopted the limited-memory Broyden–Fletcher–Goldfarb–Shanno (LBFGS) solver for weight optimization and logistic sigmoid as the activation function. The modeling procedures were implemented in Python with the open-source scikit-learn package on the web-based Google Colab platform.

Models Evaluation

The models developed based on the training dataset were also applied to the independent validation set. To assess the performances, Pearson’s correlation coefficients (r) were computed for both datasets to measure the strength of relationships between the predicted and observed values. Two metrics, including root-mean-square error (RMSE) and mean absolute error (MAE), were used to measure the deviations of the results (Equations 4 and 5). Furthermore, to obtain a scale-independent error measure in relative terms, symmetric mean absolute percentage error (SMAPE) was also computed (Equation 6) (Chen et al., 2017). The measure is resistant to outliers and has both the lower (0%) and the upper (200%) bounds.

\begin{array}{l} R M S E = \sqrt{\frac{1}{n} \sum_{i = 1}^{n} {({\hat{y}}_{i} - y_{i})}^{2}} & (4) \end{array}

\begin{array}{l} M A E = \frac{1}{n} \sum_{i = 1}^{n} | {\hat{y}}_{i} - y_{i} | & (5) \end{array}

\begin{array}{l} S M A P E = \frac{1}{n} \sum_{i = 1}^{n} \frac{2 \times | {\hat{y}}_{i} - y_{i} |}{| y_{i} | + | {\hat{y}}_{i} |} \times 100 % & (6) \end{array}

Where ŷ_i is the predicted value and y_i is the observed value of data point i, n is the number of data points.

Finally, the developed models were applied to all Sentinel-2 images to produce spatial distribution maps of the water quality parameters within the entire temporal range. The resulting maps were exported as assets in GEE. Using functions available in GEE, time-series charts and maps were generated to visualize the spatial and temporal patterns of each water quality parameter. To further exploit the potential of GEE, an online application was created through Earth Engine Apps (Kumar and Mutanga, 2018) which allowed users to explore and analyze water quality trends with a simple web interface. It should be noted that since ANN is not natively available in GEE, all the weights and biases of each neuron in the trained ANN models were manually transferred using GEE Python API.

Figure 2 shows the methodological flow of this study. It is noteworthy that while the present study focused on the end-to-end framework, the proposed workflow consisted of many interrelated components, including but not limited to the pre-processing steps and regression algorithms, which would simultaneously influence the modeling results and require further investigations.

FIGURE 2

Figure 2 Flowchart of the proposed framework to estimate water quality parameters using Google Earth Engine (GEE) cloud computing platform and artificial neural network (ANN) algorithm. (MNDWI, Modified Normalized Difference Water Index).

Results and Discussion

Evaluation of ANN Regression Models for All Water Quality Parameters

Table 4 shows the evaluation metrics for the ANN models, including both model calibration and validation sets, developed for each of the 22 water quality parameters. The ANN models, optimized through cross-validations, selected a wide range of numbers of input variables and neurons for different parameters. For instance, the models developed for salinity and turbidity used 8 input variables and 2 neurons in a single hidden layer, indicating that the dependent variables correlated well with the predictors and can be estimated using relatively simple networks. Contrastingly, the modeling of parameters such as silica and Chl-a utilized 9–11 variables with 6–8 neurons in two hidden layers. Despite the complex network, the modeling results were promising for these parameters. When the training dataset was considered, correlation coefficients greater than 0.8 were obtained for salinity, SS, nitrate nitrogen, total inorganic nitrogen, silica and Chl-a. Among these six parameters, five of them (except SS) could produce r ≥ 0.845 in the validation set as well. This proved the applicability of ANN to develop regression models and predict the values of particular water quality parameters in this setting.

TABLE 4

Table 4 Modeling results of the 22 water quality parameters calculated using artificial neural networks (ANN) on both training and validation datasets.

The water quality parameters investigated in this study could be categorized into physical properties, nutrient and inorganic constituents, and biological and microbiological examination. Results of the ANN models revealed that there was no clear pattern on the superiority of predicting particular groups of parameters. The physical properties of waters generally had moderate relationships with the spectral reflectance and the r values lied between 0.430 (pH) and 0.845 (salinity). The group of nutrient and inorganic constituents resulted in a wider range of r values from 0.228 (unionized ammonia) to 0.925 (silica). For biological and microbiological constituents, while Chl-a concentration led to a prediction model with promising fit (r = 0.898), E. coli and Faecal Coliforms generated the poorest results among all parameters tested in this study, with very weak correlations (r = 0.130–0.138) and large deviations (SMAPE = 154.8–157.8%) from the actual values. Although some scholars had reported substantial performance when applying regression models to estimate non-optical parameters including dissolved oxygen (Pizani et al., 2020), the agreement was not found in this study (r = 0.440). This could be explained by the complexity of coastal waters and the lack of consideration of interrelated optical properties among parameters (Kim et al., 2020).

From remote sensing perspective, some water quality parameters were considered as optically active components (Sobel et al., 2020), such as Chl-a and SS, which cause changes to the spectral properties of reflected light and are thus remotely detectable. Reviews by Gholizadeh et al. (2016) and Wang and Yang (2019) suggested that most of the existing studies concentrated on these optically active variables while other parameters without explicit optical properties were not well investigated. This study attempted to apply the same regression approach to all parameters and revealed a significant difference in the model performances in terms of correlation and other evaluation metrics. The findings in this study did not support that water quality parameters related to physical properties of water would produce better prediction performances compared to those related to nutrient and inorganic constituents. Instead, the potential to predict the levels of several types of nutrients and chemicals such as total nitrogen (r = 0.768), orthophosphate phosphorus (r = 0.712) and silica (r = 0.925) was also suggested. Nevertheless, since these water quality parameters were seldom studied using remote sensing methods (Gholizadeh et al., 2016), their spectral responses and modeling potentials would need to be further investigated.

Based on the ANN modeling results, the following sections would focus on three water quality parameters, namely Chl-a, SS and turbidity, which had relatively accurate estimations and strong relationships with the predictor variables. This also allowed in-depth analysis of the remote sensing-based techniques since these parameters have been investigated in similar studies.

Variables Used in Chl-a, SS and Turbidity Models

The variables used as inputs in the ANN models for Chl-a, SS and turbidity are described in Table 5. Through cross-validations, each model selected 8–11 predictors which provided relatively high correlations with the dependent variables. For the Chl-a model, it adopted a combination of band 1 (coastal aerosol), band 2 (blue), band 8A (NIR), blue-green ratio, three-band ratio of band 2, band 3 (green) and band 4 (red), as well as line-heights of band 2, band 4 and band 8A.

TABLE 5

Table 5 The variables selected as inputs in the ANN models developed in this study for chlorophyll-a (Chl-a), suspended solids (SS) and turbidity.

As suggested in existing literature, there is little consistency in the effectiveness of different bands or band ratios in estimating various water quality parameters, which further depends on the parameter ranges and study areas (Matthews, 2011; Topp et al., 2020). For instance, although band 5 in Sentinel-2 corresponds to a well-known reflectance peak of Chl-a near 705 nm in high-biomass waters (Gons, 1999), none of the band combinations directly associated with this wavelength was selected in the model. While the three-band ratio variables were originally constructed using the red-edge wavelengths from 665 to 740 nm and were shown to be useful in Ansper and Alikas (2019) and Free et al. (2020), it became one of the significant predictors of Chl-a when the same equation was applied to the visible bands in the present study. These were coherent with previous studies of Hong Kong water, which proved that a locally calibrated model outperformed models established in other study areas including the standard Case-2 Regional/Coast Color (C2RCC) processing chain model (Hafeez et al., 2019).

Despite the variety of predictor combinations, the selected variables were consistent with the prominent spectral characteristics such as absorption and scattering features of the optical constituents. As shown in the Chl-a model, 6 out of the 9 selected variables were related to coastal aerosol, blue and green bands and blue-green ratio, which represented the theoretical behavior of Chl-a that absorbs blue lights and reflects green lights (Soriano-González et al., 2019). The properties were more evident for dinoflagellates and diatoms (Sadeghi et al., 2012), the two major phytoplankton classes present in Hong Kong waters (Cheung et al., 2021). The incorporation of band 8A also matched the observation by other studies that NIR bands could be helpful to separate Chl-a from other constituents in turbid coastal waters (Hafeez et al., 2019; Ouma et al., 2020).

For the SS and turbidity models, the selected variables were dominated by those generated from green, red and red-edge bands, their square and cubic values, as well as their line-heights against neighboring bands. Similarly, refractive index and grain size of SS in Hong Kong waters were found to be sensitive to the green and red wavelengths in previous work (Nazeer and Nichol, 2015). Red-edge wavelength also becomes influential when inorganic SS concentrations increase in turbid waters (Topp et al., 2020). For turbidity, algorithms involving the red band are generally adopted owing to the scattering of particulates (Matthews, 2011). In this study, it was interesting that 7 out of 8–11 variables appeared in both SS and turbidity models, and they used small numbers of neurons (2–5) in the single hidden layer. This provided another piece of evidence that these two water quality parameters are closely related to each other and have relatively simple relationships with spectral variables.

Performances of Chl-a, SS and Turbidity Models

Figure 3 compares the predicted and observed Chl-a, SS and turbidity values using both the training and validation datasets. While the ANN models performed well in the training phase as shown in the correspondence between the fitted lines and the one-to-one lines, overfitting was also not a problem in the models since they generated similar performances when applied to the independent validation sets. Considering the validation performance, for Chl-a, the ANN model was able to accurately predict a data point with high observed value, which likely represents an occasional algal bloom event. Further analysis using Spearman’s correlation showed that the correspondence between the predicted and observed values was also satisfactory without the effect of this outlier point. However, the model tended to slightly underestimate the Chl-a concentrations when the observed values were below 4 μg/L. For SS, while the fitted line corresponds well with the one-to-one line, Pearson’s correlation coefficient was lower (r = 0.653), probably due to the larger spread of values around the fitted line (RMSE = 4.513 mg/L, SMAPE = 76.4%). In contrast to Chl-a, the model tended to overestimate the SS values when the actual values were low (< 4 mg/L) and this pattern was visible in both the training and validation datasets. For turbidity, the agreement was also strong in terms of correlation coefficient (r = 0.701) and other error metrics (RMSE = 1.954 NTU, SMAPE = 40.3%). Nonetheless, when the observed values were large, overestimation could be an issue found from the scatterplot.

FIGURE 3

Figure 3 Scatterplots comparing the predicted and observed values of chlorophyll-a, suspended solids and turbidity using the training set (A–C) and validation set (D–F). The gray dotted lines indicate the one-to-one relationships and the red lines indicate the best-fitting lines. Note that the x- and y-axes are in logarithmic scale, while the best-fitting lines and all statistics are calculated based on the untransformed scale. Besides Pearson’s correlation coefficient (r) and root-mean-square error (RMSE), Spearman’s rank correlation (ρ) is also computed to provide a measure with lower sensitivity to outliers in the dataset.

The prediction errors for Chl-a found in this study (55.1%) exceeded the desired relative error of 35% set by NASA’s Ocean Biology and Biogeochemistry Program (McCain et al., 2006). However, the standard was specified for open ocean dominated by phytoplankton. Coastal and estuarine waters are commonly referred to as Case-2 waters (Morel and Prieur, 1977), with optical properties significantly influenced by various constituents including yellow substance and suspended materials, thus were expected to be challenging to develop algorithms for predicting water quality parameters.

To serve as a benchmark of the performance, modeling results were extracted from previous works on Hong Kong water and the accuracies were found to vary for different parameters (Table 6). For Chl-a, the results were comparable and had similar levels of correlations and errors. While the prediction of SS seemed to be less satisfactory than those obtained from previous methods, the prediction of turbidity produced lower correlation but lower errors as well. This comparison indicated that, in addition to the modeling algorithm, the selection of input data and pre-processing steps were more influential on the results.

TABLE 6

Table 6 Comparison of modeling results with other studies focusing on Hong Kong waters.

In view of this, the findings presented in this study were encouraging considering the following aspects: (i) This study represented an operational context that extracted all available images acquired by Sentinel-2A/B satellites, matched with station-based measurement and applied a fully automatic analysis workflow with the capability of GEE; (ii) While linear regression analysis is more popular in GEE due to easy implementation (Tamiminia et al., 2020), the use of ANN algorithm through this platform enabled manipulation of non-linear relationships (Forkuor et al., 2017) and rapid estimation of water quality parameters when new images become available; (iii) This study was more robust than previous works in the same study area, due to the larger number of observations and extensive use of images from different years and months, facilitating the applicability over the spatially and temporally variable water conditions; (iv) Besides the commonly used RMSE and MAE, this study also reported relative errors which enable meaningful comparisons across models developed in different water regions (Shcherbakov et al., 2013).

Visualization of Model Performance Using Time-Series Charts

The developed models were applied to all Sentinel-2 images to produce spatial distribution maps of Chl-a, SS and turbidity from 2015 to 2021. To visualize the temporal patterns of these water quality parameters, time-series charts were generated by plotting the average values estimated from each remotely-sensed image against time (Figure 4). The line graph illustrated that the shapes and magnitudes of remote sensing estimates generally matched well with station-based measurements. For example, the temporal changes of Chl-a concentration showed higher values commonly found during summer, which were revealed from both remote sensing estimation and station-based measurement. From the line graph, some of the peaks found from station-based measurements could also be obtained using remote sensing estimation, showing that the method proposed in this study could capture some of the extreme events as well. It should also be noted that differences in sampling locations and frequencies of the two measurement methods could affect their correspondence presented in the graph.

FIGURE 4

Figure 4 Time-series line chart of (A) chlorophyll-a concentration, (B) suspended solids, and (C) turbidity in Hong Kong water from 2015 to 2021, obtained from station-based measurement (blue lines) and remote sensing estimation (orange lines) respectively. For station-based measurement, values are calculated as the mean value of all measurements in all locations taken in the same month. For remote sensing estimations, values are calculated as the mean of all estimated values at the locations of monitoring stations in each image. Note that the two measurement methods also differ in sampling locations (e.g. affected by cloud cover and water pixels) and frequencies (e.g. image acquisition dates).

As shown in Figure 4, the temporal range in this study was separated into model development (2015–2019), validation (2020) and prediction (2021) according to the developed methodology. Historical observations are valuable in time-series to estimate the model coefficients, followed by predicting new values and identifying abnormal changes in the monitoring period (Hamunyela et al., 2020; Elhag et al., 2021). This was also demonstrated in this study and water qualities in later years were estimated using models calibrated from previous data. Since the data release of station-based measurement by the Hong Kong government usually has a time lag of more than half a year, continuous monitoring of the latest year (2021 in this study) could be complemented by the remote sensing approach.

Spatial and Temporal Variations of Chl-a, SS and Turbidity Across Hong Kong Waters

The spatial and seasonal variabilities in the whole study area are also visualized through high-resolution maps (Figure 5). In all seasons, it showed higher SS loads and turbidity in the western part compared to the eastern waters, mainly due to the contrasting effects of Pearl River discharge from the west and oceanic currents from the east (Lai et al., 2016). High values particularly during spring and summer, as a consequence of the southwest monsoon and strong riverine inputs (Environmental Protection Department, 2021), were also evidenced in the maps. The high Chl-a concentrations in the shallow waters of Deep Bay were consistent with those revealed in Hafeez et al. (2019) (8.0–15.6 μg/L), owing to domestic wastewater and agricultural runoff (Zhou et al., 2007).

FIGURE 5

Figure 5 Spatial distribution maps of chlorophyll-a concentration, suspended solids and turbidity in Hong Kong water in the four seasons, calculated as the mean value of all remote sensing-based estimations within each meteorological season in 2015–2021. (DB, Deep Bay; VH, Victoria Harbour; TH, Tolo Harbour; SW, Southern Water).

In Victoria Harbour, while moderate values of Chl-a (7.0–10.0 μg/L) were reported in Hafeez et al. (2019), this study tended to observe lower concentrations on average (2.0–4.0 μg/L). This was probably because the present study utilized a larger number of images and Chl-a concentration tends to remain at low levels with high flushing rates most of the time throughout the year (Xu et al., 2010). Tolo Harbour, a land-locked inlet located in northeastern Hong Kong, was suggested to have highly eutrophic water and regular algal blooms due to restricted water exchange, low tidal flushing and persistent stratification (Lai et al., 2016). The maps produced in this study provided further observations that high Chl-a concentrations in that area tended to occur more frequently in autumn and winter. This seasonal variation matched with the findings of a recent forecasting study by Deng et al. (2022) and could be related to the downwelling process and concentration of algal biomass in that period (Xu et al., 2010). In the southern waters, slightly higher Chl-a concentrations were observed in summer from the produced maps, which also agrees with the observations by Zhou et al. (2012) that these areas represent the estuarine coastal frontal region and favor the development of algal blooms in July and August.

Besides seasonal trends, the robustness of the developed model was further demonstrated when applied on a particular occasion. During a red tide incident that occurred in Victoria Harbour from 13^th to 16^th July 2020 caused by phytoplankton species Amphora spp. (Agriculture, Fisheries and Conservation Department, 2021), the Chl-a concentration map produced from the Sentinel-2 image acquired on 13^th Jul 2020 successfully revealed high values in that location (Figure 6). The estimated Chl-a values in the nearby areas were also consistent with station-based measurements within the same period. While moderate to high concentrations were found in Victoria Harbour and North Western Water Control Zones, the map was able to present low concentrations, on the contrary, in Port Shelter at the same period. These results again confirmed the prediction ability in the monitoring period using models calibrated from historical data, where potential anomalies in water conditions including algal bloom events and their spatial variabilities could be captured.

FIGURE 6

Figure 6 (A) True color composite of Sentinel-2 image acquired on 13^th July 2020 and (B) the estimated chlorophyll-a (Chl-a) concentration map, overlaid with in situ data measured in the monitoring stations.

Customized GEE Application

Finally, using functions available in GEE, an online application was created (https://khoyinivan.users.earthengine.app/view/marine-water-quality-hk) to allow users to explore and analyze water quality trends with a simple web interface (Figure 7). Besides zooming and panning the map to show different locations, users can interact with the applications in several ways, including displaying the concentration map of water quality parameters for a specific date and plotting time-series charts of parameters for a specific location.

FIGURE 7

Figure 7 A screenshot of the developed online application in Google Earth Engine showing the water quality map for a specific date and time-series charts for a specific location.

Achievements and Limitations of This Study

The online application, together with the analysis workflow, takes advantage of the cloud computing GEE infrastructure and the readily available Sentinel-2 satellite images. In this study, the comprehensive archive of Sentinel-2 dataset was seamlessly searched and processed in the cloud platform, without the need for downloading image data as in the conventional processing approach (Kumar and Mutanga, 2018; Li et al., 2021). The method is cost-effective as the Sentinel-2 data, GEE and Python software are open-source and freely available. Sentinel-2 will revisit the study area every five days with the two-satellite constellation and become available in GEE (Rudiyanto et al., 2019). Since the mapping processes can also be implemented automatically with the algorithm coded in GEE and Python, this enables monitoring from continuously updated remote sensing datasets in a timely manner and generates actionable information for decision-making (Markert et al., 2018; Page et al., 2019). Although the obtained levels of accuracies were below the desired standard, the end-to-end cloud computing workflow demonstrated in this study should be investigated considering the cost and computational efficiency for rapid information delivery.

The developed methodology and GEE scripts can further be extended to larger scales to increase usability and transferability (Kumar and Mutanga, 2018). It was anticipated that integrated use of remote sensing, in situ measurements and computer modeling can lead to an enhanced knowledge of the water quality (Gholizadeh et al., 2016). Still, the remote sensing tool would require support from existing in situ monitoring infrastructure and established laboratories (Pompêo et al., 2021). Hence, the routine measurements of biophysical and chemical properties of water by EPD of the Hong Kong government provide great opportunities for experimentation in this study.

Nevertheless, the remote sensing approach is limited by satellite repeat cycles and cloud coverage. Sent et al. (2021) described that in their work, only 55% of the total images available over the 2-year study period were cloud-free and in good atmospheric conditions. The percentage was even smaller (31.7%, 120 out of 379 images) in this study (Table 3), due to the subtropical oceanic monsoon climate in Hong Kong (Tian et al., 2014). Moreover, the image availability was non-homogeneous throughout the year with fewer images available in the wet season (spring and summer), making routine water monitoring even more challenging. Future work can explore the combined use with other satellite missions, such as Landsat satellites (Page et al., 2019), to enhance the observing capability.

Another limitation in the model development process in this study was the lack of in situ reflectance data. The calibration and evaluation of the algorithms were performed solely using satellite-based reflectance and station-based water quality measurements, which inhibited the understanding of the atmospheric influence and the true relationships between water-leaving reflectance and water quality parameters (Yoon et al., 2019). As illustrated in Kuhn et al. (2019), 70–90% of the signals in water pixels came from atmospheric effects and the choice of correction methods could bias parameter retrieval accuracies, posing challenges to land observation satellites like Sentinel-2. Adjacency to land areas further affected the optical complexity of coastal waters and hindered the use of algorithms developed for open oceans (Moses et al., 2017), as evidenced by the high NIR and SWIR reflectance of water pixels in this study. Many scholars have performed intercomparison of various atmospheric correction and bio-optical algorithms (Gernez et al., 2017; Poddar et al., 2019; Soomets et al., 2020; Sent et al., 2021), using both in situ and satellite measurements, to assess the retrieval performances in specific regional waters. Future work can focus on the identification of optimal algorithms for predicting water quality parameters, including both optically active and inactive parameters, that can be incorporated into GEE applications for operational monitoring.

Conclusion

Using coastal waters of Hong Kong as the study area, this study utilized GEE to (i) query and pre-process all Sentinel-2 observations that coincided with in situ measurements; (ii) extract the spectra to develop empirical models for different water quality parameters using ANN; and (iii) present the results using spatial distribution maps, time-series charts and an online application. Together these features enabled novice users to discover, analyze and visualize spatial and temporal patterns from remotely-sensed images. This study demonstrated that the integration of Sentinel-2 satellite images, in situ measurements and the GEE cloud computing platform can offer new opportunities for providing cost-effective and timely prediction of water quality parameters. With multiple water quality parameters analyzed in this study, the mapping results could help comprehensive understanding the dynamics of water qualities within the study area in multiple years, raise awareness of environmental changes in relation to human activities and promote sound management policies. When longer records of both water quality measurements and freely available earth observation data become available in the future, more integrated and robust predictive models can be developed to continue to contribute towards an effective water resource monitoring framework.

Data Availability Statement

Publicly available datasets were analyzed in this study. This data can be found here: https://github.com/ivanhykwong/marine-water-quality-time-series-hk.

Author Contributions

IK and FW initiated this study. IK designed the methods, conducted the analysis, and drafted the manuscript. FW and TF supervised the work and improved on the result interpretation. All authors contributed to the article and approved the submitted version.

Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s Note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Acknowledgments

The authors would like to thank the Environmental Protection Department (EPD) of the Hong Kong government for providing station-based coastal water quality data. Sentinel-2 image archive was provided by the European Space Agency (ESA) and acquired through the Google Earth Engine API. We thank the editors and all reviewers for their valuable comments on this manuscript.

References

Agriculture, Fisheries and Conservation Department (2021) Hong Kong Red Tide Database. Available at: https://redtide.afcd.gov.hk/urtin/#/maps?mode=1&lang=en (Accessed 31 January 2022).

Google Scholar

Amani M., Kakooei M., Moghimi A., Ghorbanian A., Ranjgar B., Mahdavi S., et al. (2020). Application of Google Earth Engine Cloud Computing Platform, Sentinel Imagery, and Neural Networks for Crop Mapping in Canada. Remote Sens. 12 (21), 3561. doi: 10.3390/rs12213561