A Multivariate Model for Coastal Water Quality Mapping Using Satellite Remote Sensing Images.

his study demonstrates the feasibility of coastal water quality mapping using satellite remote sensing images. Water quality sampling campaigns were conducted over a coastal area in northern Taiwan for measurements of three water quality variables including Secchi disk depth, turbidity, and total suspended solids. SPOT satellite images nearly concurrent with the water quality sampling campaigns were also acquired. A spectral reflectance estimation scheme proposed in this study was applied to SPOT multispectral images for estimation of the sea surface reflectance. Two models, univariate and multivariate, for water quality estimation using the sea surface reflectance derived from SPOT images were established. The multivariate model takes into consideration the wavelength-dependent combined effect of individual seawater constituents on the sea surface reflectance and is superior over the univariate model. Finally, quantitative coastal water quality mapping was accomplished by substituting the pixel-specific spectral reflectance into the multivariate water quality estimation model.

Introduction

In light of the synoptic spatial coverage and routine availability of satellite images, there have been many implementations of environmental monitoring using remote sensing images in recent years. Among these practices, applications focusing on terrestrial and water environments such as landuse/landcover classification and change detection, landslide site detection, lake and reservoir trophic state monitoring [1-8], ocean harmful algal bloom monitoring [9-12], and coastal water quality monitoring [13-21] are most frequently conducted, although there have been significantly fewer examples of coastal water quality monitoring, as compared to monitoring of impounded water bodies such as lakes and reservoirs. In our opinion, it may be attributed to several factors including difficulty in field data collection, less human activities in the target coastal areas, and the public attitude towards coast utilization and management.

An essential difference in the nature of water quality in coastal area and impounded water bodies is that water quality of impounded water bodies is almost never affected by the downstream condition, whereas coastal water quality is constantly affected by upstream flow discharge and seawater flushing. Upstream inflow may contain high concentration sediments, especially during high-flow periods. Suspended solid material trapped in impounded water bodies can only be removed through overflow from the outfall channel or by other means of surface discharge. For lakes of no flow outlet, suspended solid material carried in by upstream inflow may have long-lasting effect on lake water quality. On the contrary, suspended solid material in coastal water may be mixed or flushed by seawater, and thus the characteristics of creeping impact and low immediate effect make it difficult to sense the emerging consequences which may be severely deleterious and irreversible on coastal ecosystems. For example, Tomascik and Sander [22] and Hoegh-Guldberg et al. [23] reported changes to coastal discharge as one of the most serious threats facing coral reef ecosystems worldwide. Thus, understanding the long term effect of sediments carried in the upstream discharge on coastal water quality necessitates a routine monitoring scheme.

Efficient monitoring of the coastal water quality cannot be achieved using only water samples taken in vast coastal area. Not only is the water sampling time and labor consuming, but the required number of samples would also practically hinder such implementation. An integration of remote sensing technique and seawater sampling is therefore more appealing and worth pursuing. Thus, the purpose of this study is to investigate the feasibility of using satellite remote sensing images for coastal water quality monitoring. Specifically, generation of distribution maps of three water quality variables including the total suspended solid (TSS), turbidity (Tb), and Secchi disk depth (SDD) from satellite images is desired.

Study area and data collection

A coastal area near the outlet of a flood diversion tunnel (the Yuan-Shan-Tzu Diversion Tunnel) in northern Taiwan was chosen for this study (Figure 1). The YST tunnel was constructed to divert flood flow from the upper Keelung River Basin to a discharge outlet at the northern tip of the island. The YST tunnel, completed in 2003 with a diameter of 12 m and 2.48 km in total length, is capable of diverting approximately 81% (1,310 m3/s) of the 200-year flood flow (1,620 m3/s) at a cross section near the inlet of YST tunnel. Within a five month period (July to November of 2007), a few water sampling campaigns were conducted in a coastal area north of the tunnel outlet. The sampling dates and relevant storm information are shown in Table 1. Water samples were taken within 0 – 20 cm range below the sea surface at eight locations (see Figure 2) during each sampling campaign. Global positioning systems were used to guide the sampling vessel to the desired sampling locations. Considering the west-to-east surface current direction of the season, the sampling area extends from a little northwest of the outlet to about 2 km to the east of the outlet. Sampling point 4 is located within an area known as the Yin-Yang Sea. Geology of the near Yin-Yang Sea area has a large amount of pyrite that does no dissolve easily in water. The Yin-Yang Sea area, just offshore from an old metal mining township, frequently receives runoff containing high iron ion concentration, making the sea surface visually distinct. Secchi disk depth and turbidity were measured in situ and water samples were taken to the Environmental Chemistry Lab at the National Taiwan University for analyses of suspended solids and other properties. Table 2 summarizes statistical properties of the three water quality variables. Two sampling campaigns (09/20/2007 and 10/08/2007) took place one day after activation of flood diversion. Comparison of the water quality data of the no-diversion and post-diversion periods is shown in Figure 3. Differences in medians and ranges of TSS, Tb and SDD are apparent. For example, median of SDD drops from 6.8 m of the no-diversion period to 3.8 m of the post-diversion period, whereas median of TSS increases from 1.6 mg/L of the no-diversion period to 5.2 mg/L of the post-diversion period. Also, excluding the outliers, the ranges of the water quality variables of the no-diversion and post-diversion periods are almost non-overlapping.

Figure 1.

Location map of the study area. Dashed circle indicates the coastal area where water sampling was conducted. Shaded area shows the Keelung River watershed.

Table 1.

Dates of water sampling and SPOT image acquisition.

Sampling date	SPOT image acquisition date	Relevant storm events	Volume of diverted flow (m3)
7/02/2007	7/04/2007(SPOT-4)	No storm	0
7/18/2007	7/19/2007(SPOT-4)	No storm	0
8/15/2007	NAa	No storm	0
8/23/2007	8/23/2007(SPOT-5)	Typhoon Sepat(8/16∼8/19)	0
9/07/2007	9/03/2007(SPOT-5)	No storm	0
9/20/2007	NAa	Typhoon Wiphab(9/17∼9/19)	1,051,200
10/08/2007	NA	Typhoon Krosab(10/4∼10/7)	16,133,400
11/14/2007	NA	No storm	0

Satellite images were not collected due to high percentage of cloud cover.

Flow diversion activated.

Figure 2.

Water sampling locations in the coastal area near the YST tunnel outlet (a) and the radiometric control area (b).

Table 2.

Statistical properties of water quality variables.

	Mean	Standard deviation	Maximum	Minimum
Secchi disk depth (m)	5.39	2.11	10.30	0.50
Turbidity (NTU)	2.19	4.19	29.50	0.38
Total suspended solid (mg/L)	4.79	5.95	30.61	0

Total number of samples: 62

Figure 3.

Box plots of water quality data of no-diversion and post-diversion periods.

Before pursuing establishment of water quality estimation models using the water quality data and remote sensing images, we conducted a careful check on measurements of water quality variables. The purpose of such data check is to screen out data which might have been contaminated by inappropriate sampling of water samples or erroneous measurement in the lab. In general, the Secchi disk depth, total suspended solids, and turbidity are inter-related, as demonstrated in Figure 4. The points marked by dashed-circles are significantly inconsistent with such correlation, and thus are excluded in subsequent analyses. Koponen et al. [24] also showed similar relationship between the Secchi disk depth and the turbidity.

Figure 4.

Empirical relationships among different water quality parameters. Points marked by dashed-circles are not included in regression modeling and are excluded from the subsequent analyses.

A few multispectral images from SPOT satellites with acquisition dates close to the dates of sampling campaign were also collected (see Table 1). During and immediately after the Wipha and Krosa typhoon events, the study area were almost completely under cloud cover, and thus no SPOT images were collected. The multispectral SPOT images include images of three spectral bands – green (0.5–0.59μm), red (0.61–0.68μm), and near infrared (0.78–0.89μm), with pixel resolution of 20 m for SPOT-4 or 10 m for SPOT-5.

Remote sensing image analysis

The amount of solar radiance reaching the satellite sensor can be expressed as [25]: where

L= solar radiance reaching the sensor

= the exoatmospheric solar irradiance

E = the downwelled irradiance from the sky dome onto the target

L= upwelled solar radiance

θ = the view angle in the sensor-target direction

σ = the sun angle in the sun-target direction

ϕ = the azimuth angle in the sun-target direction

ϕ = the azimuth angle in the sensor-target direction

τ1(λ) = the atmospheric transmittance along the sun-target path

τ2(λ) = atmospheric transmittance along the target-sensor path.

λ = spectral wavelength of solar radiation

ρ = spectral reflectance of the target surface

F = the obstruction factor.

The obstruction factor F in Equation (1) accounts for the proportion of irradiance that may be obstructed by adjacent objects or surface slope of the target. If the target object is on a horizontal surface and free of adjacent object obstruction, the factor F equals 1. It is also worthy to note that the sun and view angles are defined with reference to the normal of the target surface. If the target is located on a slope, the sun and view angles will need to be adjusted accordingly. Readers are referred to Schott [25] for detailed calculation of solar radiances arriving at the sensor.

The reflectance ρ varies with spectral wavelength and orientation angles. If the target object is assumed to be a diffuse reflector with a constant reflectance ρ (λ) in all directions, we then have

On the right-hand side of the above equation, only the spectral reflectance ρ(λ) represents physical property of the target surface. The upwelled radiance L does not even get into contact with the target.

Environmental monitoring using remote sensing images often requires derivation of physical properties (reflectance, for example) of the target objects from satellite images. Unfortunately, the upwelled radiance L, the atmospheric transmittance τ1 and τ2, the downwelled irradiance E, and the exoatmospheric solar irradiance are generally not available for most applications, and we have to resort to other means for estimation of the reflectance.

For most local-scale environmental monitoring applications, L, τ1, τ2, E, and can be assumed constant (or spatially invariant) within the study area. While on the contrary, the sun angle σ and the obstruction factor F are dependent on the surface slope of the target, and the reflectance ρ(λ) is dependent on surface cover of the earth. Their values may vary from pixel to pixel within a scene. If only pixels on horizontal surface and free of adjacent obstruction are considered (F = 1), Equation (1) may be expressed as: where and k2 = L(θ, ϕ).

A common practice dealing with the upwelled radiance L(θ, ϕ) in satellite remote sensing is the dark object subtraction (DOS) method [5, 26, 27]. The basic concept of the DOS method is to identify very dark features within the scene. The minimum scene radiance is set to be the upwelled radiance based on the assumption that it represents the radiance from a pixel with near zero reflectance. If the minimum scene radiance is subtracted from the radiance of each individual pixel, the processed image is then assumed free of atmospheric scattering effect.

After removing the upwelled radiance L using the DOS method, the DOS-adjusted radiance is linearly related to the surface reflectance, i.e.:

Based on the linear relationship between and ρ(λ), we devise a surface reflectance estimation scheme through reflectance calibration in a radiometric control area (RCA).

In this study a radiometric control area of approximately 30 m × 60 m was chosen for spectral reflectance calibration. The RCA is a horizontal paved open area with homogeneous and stationary surface reflectance and no adjacent obstruction (see Figure 2). It is located in a restricted and free of public access harbor area. The wavelength-dependent surface reflectance ρ(λ) of RCA is then calibrated using a variable spectral radiometer (VSR) which is equipped with two spectral-variable filters capable of detecting spectral radiances in various 7 nm-wide windows within the 0.4 – 0.72 μm and 0.65 – 1.1 μm ranges. The VSR was moved around within the radiometric control area taking multispectral images. When taking images within the RCA, a standard reflectance disk which has been pre-calibrated to have ρ(λ) ≈ 1 over the 0.25 – 1.1 μm wavelength range was also placed within the viewing area. Reflectance of the radiometric control area is then calculated as the ratio of average radiance from RCA to average radiance from the standard reflectance disk, i.e.: where ρ(λ) is the reflectance of RCA, and L̄ and L̄ are respectively average radiances received at VSR sensor from the RCA surface and from the standard reflectance disk. For RCA reflectance calibration, the effect of upwelled radiance can be neglected since the VSR is placed near the ground surface. The RCA-average reflectances with respect to various spectral wavelengths are shown in Figure 5. The RCA-average reflectances corresponding to the green, red, and near infrared SPOT spectral bands (hereafter referred to as the RCA band reflectances) are calculated to be 0.097, 0.113, and 0.161%, respectively. The RCA band reflectances are considered constant since the land surface condition within the RCA is relatively homogeneous and stationary.

Figure 5.

Calibrated RCA-average reflectances and band-average reflectances.

Assuming the sea surface is horizontal, the DOS-adjusted radiances of pixel A on the sea surface and pixel B within RCA are respectively expressed by: and

Thus, reflectance of the sea surface pixel A can be calculated as

In practice, and ρ(λ) are respectively substituted by the average radiance and reflectance of the RCA, i.e. where represents the average value of DOS-adjusted radiances within the RCA and ρ(λ) is the RCA band reflectance. The reflectance calibration ratio in the above equation may vary among different SPOT scenes since varies due, to scene-to-scene variations in orientation angles and atmospheric transmittance. Table 3 summarizes the reflectance calibration ratios of individual SPOT multispectral images.

Table 3.

Scene reflectance calibration ratios of SPOT multispectral images used in this study.

Image acquisition date	Reflectance calibration ratio ρ_c(λ)L¯sλ,c′(θ,ϕ_d)
	Green	Red	Near infrared
7/04/2007	0.00282	0.00313	0.00420
7/19/2007	0.00232	0.00249	0.00410
8/23/2007	0.00331	0.00336	0.00532
9/03/2007	0.00345	0.00342	0.00500

Water quality estimation

In order to map the spatial distribution of the water quality variables using remote sensing images, it is necessary to establish water quality estimation models based on the reflectance of sea surface. A few simple or multiple regression models have been used in the literature [5-8, 13, 14, 21, 28-37]. Most of these models fall into one of the following forms: where Y represents a water quality variable and X can be reflectance of a specific spectral band, ratio of reflectances of different spectral bands, or other arithmatic calculation of band reflectances.

In order to choose appropriate models for water quality mapping, we first examined the scatter plots of the water quality measurements versus the band-dependent sea surface reflectances, as shown in Figures 6(a)-(c). Although the data points are widely dispersed, particularly in the lower measurement ranges, measurements of the turbidity and the total suspended solids tend to increase with the sea surface reflectance, whereas the secchi disk depth tends to decrease with increase of the sea surface reflectance. To better illustrate the variation trend of these water quality variables, the water quality measurements were grouped into several incremental intervals. The interval-average water quality and the corresponding average sea surface reflectance were then calculated. Figure 6(d) demonstrates that the interval-average water quality variables are well related to the red band reflectance. Using Equation 10(a), the specific water quality estimation models are as follows: where SDD, Tb, and TSS are respectively measured in units of m, NTU, and mg/L, and ρ is the red band spectral reflectance in percentage. The signs of regression coefficients of ρ in Eqs. (11)–(13) are consistent with the physical phenomena normally observed in the natural environment.

Figure 6.

Scatter plot of water quality measurents versus band-dependent sea surface reflectances. (a), (b), (c) Original water quality measurements, (d) interval-average water quality measurements.

Although the above water quality estimation models were established using the interval-average measurements, these models were adopted for water quality estimation using the pixel-based sea surface reflectance. The applicability of these models was checked by comparing the model estimates against the original measurements. Figure 7(a) demonstrates that the water quality estimates are roughly consistent with the corresponding measurements. Each of the above models utilizes the single band reflectance ρ for estimation of a single water quality variable, and is referred to as the univariate model in this study.

Figure 7.

Water quality measurements versus estimates using the univariate and multivariate estimation models. [Note that there is an out-of-bound SDD estimate by the univariate model.]

The water body is a mixture of seawater, suspended solids, dissolved organic matters, zooplankton, etc. The sea surface reflectance of a specific wavelength is affected by the combined effect of these constituents. On the other hand, the effects of individual constituents on the sea surface reflectance vary among different spectral wavelengths. Such wavelength-dependent combined effect must be reflected in the water quality estimation model. Thus, we propose the following multivariate model for water quality estimation using the multispectral reflectances: where β's represent the regression coefficients, ρ and s are respectively the sea surface reflectance and the random component with respect to the spectral band or the water quality variable shown in their subscripts, and n is the number of water quality measurements.

The above multivariate regression model can also be expressed as where Y, X and W are respectively matrices of the water quality measurements, the surface reflectance, and the regression coefficients, and E is the error matrix. The least squared estimator of the regression coefficients is [38]:

Using a total of 25 samples of water quality measurements and their corresponding multispectral sea surface reflectances, the regression coefficient matrix W is estimated to be

Thus, water quality estimates of the multivariate model are obtained by the following equation:

Figure 7(b) demonstrates the results of water quality estimation using the multivariate model. Comparing to the estimation results of the univariate model, water quality estimates of the multivariate model show significantly less degree of dispersion around the line of equivalence. Notably, there is an out-of-bound estimate of SDD using the univariate model (see Figure 7(a)), whereas the corresponding estimate by the multivariate model is quite accurate. The superior of the multivariate model may be attributed to its capability of modeling the wavelength-dependent combined effect of the seawater constituents. It is also noteworthy that the sign of regression coefficients of the red band reflectance (ρ) in both models are consistent, and the red band reflectance remains the dominant factor in the multivariate model, judging from its significantly highest absolute value of the regression coefficient.

Water quality mapping

Quantitative coastal water quality mapping was accomplished by substituting the pixel-specific spectral reflectance calculated by Eq. (9) into the multivariate model of Equation (18). The resultant water quality distribution maps are shown in Figures 8, 9 and 10. In general, SDD increases outward from the near shore area, whereas decreasing Tb and TSS can be observed. Such spatial variation trends are particularly evident on August 23 and September 3, 2007. Spatial variations of water quality variables on July 4 and 19 appear to be more complicated and the influence of the water quality condition in the Yin-Yang Sea area is more significant. As can be seen in Figure 11 (a), on July 2, 2007 there were higher waves on the sea and the brownish water color (the lower right corner in the figure) suggests higher concentrations of TSS and Tb in seawater. In contrast, the sea condition was calm on August 23 and no brownish water color can be observed in Figure 11(b). Turbidity distribution mapping also reveals areas with negative Tb values. Although small in magnitude, these negative estimates suggest exercising extra caution for estimates of low Tb concentration.

Figure 8.

Spatial distribution of secchi disk depth.

Figure 9.

Spatial distribution of turbidity.

Figure 10.

Spatial distribution of total suspended solids.

Figure 11.

Photos of the Yin-Yang Sea area taken during water sampling campaigns.

Within the study area, the Yin-Yang Sea area has exceptionally higher TSS and Tb concentrations and lower SDD values due to non-dissolved solids routinely received from its upstream area. As for the area offshore from the YST tunnel outlet, no significant effect of the diverted flood flow on coastal water quality has been observed. This may be attributed to very few cases and short duration of flow diversion since the tunnel completion in 2003. However, a monitoring routine using satellite images is recommended for assessing the long term effect of the diverted flood flow on the coastal water quality.

Conclusions

In this paper we demonstrate the feasibility of coastal water quality mapping using remote sensing images. A few concluding remarks are drawn as follows:

A surface reflectance estimation scheme which involves choosing a radiometric control area was proposed in this study. The scheme is applicable for local-scale environmental monitoring applications.

The three water quality variables (TSS, Tb and SDD) are found to be most related to the red band surface reflectance. High values of the sea surface reflectance generally correspond to high TSS and Tb concentrations and low SDD values.

The water body is a mixture of the seawater and other constituents including the suspended solids, the dissolved organic matters, the zooplankton, etc. The proposed multivariate water quality estimation model takes into consideration the wavelength-dependent combined effect of individual constituents on the sea surface reflectance and yields more accurate water quality estimation results.

Water quality mapping using remote sensing images shows a general pattern of increasing SDD and decreasing Tb and TSS outward from the coast. Under higher wave condition, water quality in the Yin-Yang Sea area may have more significant influence on the spatial distribution of water quality in the nearby area.

Until present, no significant effect of the diverted flow on coastal water quality has been observed due to few cases of flow diversion. However, a routine operation of coastal water quality mapping utilizing satellite images is recommended for assessment of the long term effect of the diverted flow.