The Influence of River Morphology on the Remote Sensing Based Discharge Estimation: Implications for Satellite Virtual Gauge Establishment

Abstract

Monitoring of river discharge is a key process for water resources management, soil and water conservation, climate change, water cycling, flood or drought warning, agriculture and transportation, especially for the sustainable development of rivers and their surrounding ecological environment. Continuous and comprehensive discharge monitoring was usually impossible before, due to sparse gauges and gauge deactivation. Satellite remote sensing provides an advanced approach for estimating and monitoring river discharge at regional or even global scales. River morphology is generally considered to be a direct factor that affects the accuracy of remote sensing estimation, but the specific indicators and the extent to which it affects the estimation accuracy have not yet been explored, especially for medium to small rivers (width < 100 m). In this paper, six sites with hydrological gauges in the upper Heihe River Basin (HRB) of northwestern China and the Murray Darling Basin (MDB) of southeastern Australia were selected as the study cases. River discharge was estimated from Landsat imagery using the C/M method accordingly. River gradient, sinuosity, and width were obtained from Digital Elevation Model data for each site. Global Surface Water Dataset (GSWD) was also employed for indicating the dynamic status of river morphology. A series of methods were applied to analyze the influence of river morphology on estimation accuracy qualitatively and quantitatively, based on which we established inference about the theory of selecting satellite virtual gauges (SVGs). The results confirm the feasibility of the C/M method for discharge estimation, with the accuracy affected by multiple river morphological indicators. Among them, river width was found to be the most significant one. Moreover, water occurrence and water extent extracted from GSWD also have impact on the discharge estimation accuracy. Another independent river section in MDB was set as an example to demonstrate the reasonability of the established theory. It is anticipated that this study would promote the application of remote sensing for discharge estimation by providing practical guidance for establishing appropriate SVGs.

1. Introduction

River discharge, or called streamflow, is a vital component of the hydrological cycle. Therefore, accurately monitoring river discharge is crucial for a variety of hydrological studies and applications, such as water resources management, soil and water conservation, flood/drought monitoring under a changing climate [1,2]. River discharge monitoring via traditional stream gauging networks is usually the most appropriate approach of estimating discharge, which however is limited to point locations along the river in accessible and developed reaches [3]. On one hand, gauges are distributed unevenly and sparsely worldwide. One the other hand, the number of gauges is further decreasing in recent years due to various reasons, such as extreme weather, lack of maintenance, etc. [4]. This makes it even more difficult for natural rivers to be monitored continuously and comprehensively [5].

Remote sensing technology provides an alternative and advanced discharge estimation approach due to its unique advantages in obtaining river water extent, width, water level, and other river characteristics timely and stably [6,7], which greatly extends the source of the traditional hydrological data [8]. Gleason et al. [9] suggested that for rivers with limited observations, remote sensing could act as important supplementation and expansion. While for rivers with completely no observation, it would be the core data source. Therefore, river discharge monitoring with satellite remote sensing is a burgeoning field attracting increasing attentions [10]. Since the Prediction in Ungauged Basins (PUB) was proposed in 2003 [11], satellite remote sensing has been more and more widely used for river discharge estimation through establishing virtual gauges [12]. So far, a series of remote sensing-based discharge estimation methods have been developed. Gleason et al. [9] made a comprehensive review of these methods regarding the availability of gauge data (gauged, semi-gauged, regionally ungauged, and totally ungauged basins).

The rating curve method is a traditional approach that estimates river discharge by establishing a relationship between available observed discharge and one of the remotely sensed river characteristics, such as water level [13] and water area [14]. This approach is suitable for rivers that have in situ observations, but it cannot be applied to ungauged rivers. Other studies used remotely sensed data as major input to estimate discharge based on hydraulic geometry [7,15,16,17,18,19], or flow wave propagation models [20,21]. However, some key information required by these methods, such as river bathymetry and riverbed morphology, is usually unavailable. Although it may be obtained from remote sensing techniques [22,23], extensive and various remote sensing data are needed, which may limit their applicability in some rivers.

The Calibration/Measurement (C/M) method, first proposed by Brakenridge et al. [24], can be considered as an improved rating curve approach. It estimates discharge based on the close correlation between observed discharge and the ratio of a land pixel for calibration (C) and a water pixel for measurement (M). It was first tested using Advanced Microwave Scanning Radiometer for EOS (AMSRE), and then extended to MODIS imagery for estimating river discharge for different river basins [25,26,27,28]. Due to the coarse spatial resolution of AMSRE and MODIS data, these studies only focused on large rivers. Other studies, therefore, tried to extend the C/M method to medium-to-high resolution images, such as Landsat and Sentinel-2 [29,30], for estimating discharge of medium-to-small rivers (width < 100 m). These studies demonstrate the feasibility of utilizing various remote sensing data for estimating river discharge for those rivers that have limited observations due to various reasons.

As more and more remotely sensed data are becoming available, every single river is likely to be fully covered by various data. In that case, how to locate an ideal reach for establishing satellite virtual gauge (SVG) to properly turn remote sensing signals to discharge information? Existing studies suggested that river morphology is likely to be an important factor that affects discharge estimation accuracy. For example, Durand et al. [31] compared five different methods for estimating discharge in 19 selected rivers and found obvious accuracy difference due to different river morphology. Brinkerhoff et al. [32] demonstrated the important influence of river morphology on discharge estimation using a hydraulic method. Hou et al. [12] concluded that wide channels with strong temporal variations, broad floodplains and multiple braided or anastomosing channels provided the best conditions for SVG. They tested a series of different discharge estimation methods at different river reaches globally. However, as the C/M method is not a hydraulic method, would the river morphology also affect its accuracy? How much uncertainty would be induced by river morphology in the C/M modelling process? Answers to these questions are still not clear.

This study aims to investigate the influence of river morphology on estimating discharge using remote sensing. We focus on medium to small rivers (width < 100 m) through implementing C/M method with medium-resolution Landsat imagery, considering that smaller rivers are more extensively existing and ungauged. We try to investigate how river morphological indicators, including river gradient, sinuosity, and width, affect discharge estimation accuracy through statistical analysis, based on which we try to provide guidance for selecting appropriate river reaches to establish SVGs through a case demonstration. As most of the current remote sensing-based discharge estimation studies overlooked the site selection, this study is expected to help reexamine how to carry out C/M discharge estimation more effectively, especially when remote sensing data are usually available for the whole river. As the C/M method requires a time series of discharge data as the input, we select a series of reaches located near several documented gauges as test sites. Our objectives include: (1) qualitatively and quantitatively revealing the impact of river morphology on the accuracy of C/M method, and (2) providing guidance for selecting appropriate sites for establishing SVGs based on finding of the above analysis.

2. Materials and Methods

2.1. Materials

2.1.1. Study Area

As shown in Figure 1, two study areas, the Heihe River Basin (HRB) (Figure 1b), China, and the Murray Darling Basin (MDB) (Figure 1c), Australia, both in semi-arid region, were selected. For each study area, three gauges with observed hydrological records were selected. For the HRB, Gauge Qilian at the Babao River, Gauge Zhamashike at the Heihe River, and Gauge Sunan at the Liyuan River were selected. For each gauge, four river sections (called study sites hereafter) were selected nearby. They were marked as 1-1, 1-2, 1-3, 1-4, 2-1, 2-2, 2-3, 2-4, 3-1, 3-2, 3-3 and 3-4 (Figure 1d–f). For the MDB, Edwards River at Leiwah (Gauge ID: 409035), downstream of Balranald Weir (Gauge ID: 410130), and Murray River at Below Wakool Juction (Gauge ID: 414200) were selected, along with their corresponding nearby study sites marked as 4-1, 4-2, 4-3, 4-4, 5-1, 5-2, 5-3, 5-4, 6-1, 6-2, 6-3 and 6-4 (Figure 1g–i). These sites were selected where there are no tributaries to ensure that the streamflow in these sites is consistent with that recorded at the gauges. The size of site (21 × 21 pixels) was determined considering the computation efficiency of C/M modelling. Darling River at Wilcannia Main Channel (Gauge ID: 425008) (Figure 1c) was selected to demonstrate the implications for SVG establishment.

2.1.2. Data

Observed discharge in the HRB were obtained from Li et al. [29]. Observed discharge in the MDB were obtained from New South Wales Office of Water and Victorian Department of Environment and Primary Industries, after being converted from Megaliters/day (ML/d) to m³/s (1ML/d = 0.011574 m³/s). Their hydrographs are displayed in Figure 2. The modelling period and validation period were listed in

Table 1. Modelling period and validation period of each gauge, and the corresponding amount of Landsat images.

Gauge	Modeling Period	Amount	Validation Period	Amount
Qilian	7 January 2000–22 February 2011	126	26 March 2011–2 December 2015	31
Zhamashike	7 January 2000–26 March 2011	127	11 April 2011–2 December 2015	32
Sunan	25 January 2010–12 January 2017	40	13 February 2017–30 December 2017	11
409035	14 January 2000–20 July 2007 & 18 February 2010–26 April 2011	97	21 August 2007–16 December 2009 & 1 September 2011–4 November 2011	24
410130	14 January 2000–20 July 2007 & 18 February 2010–26 April 2011	97	21 August 2007–16 December 2009 & 1 September 2011–4 November 2011	24
414200	14 January 2000–20 July 2007 & 18 February 2010–26 April 2011	97	21 August 2007–16 December 2009 & 1 September 2011–4 November 2011	24
425008	5 January 2000–26 September 2006 & 24 January 2010–26 December 2010	80	28 August 2007–20 October 2009 & 27 January 2011–11 November 2011	19

. As the MDB experienced a 10-yr dry period with very low streamflow observed (2000–2009), which was then reversed later (2010–2011) with high streamflow [33], we selected two periods for modelling and two periods for validation for the MDB sites to ensure that both high flow and low flow periods have been included for modelling and validation. For each site, we used Level 2 atmospherically corrected surface reflectance of Landsat data, which were acquired from Google Earth Engine (GEE) [34]. We checked the data for each site, and manually picked out those cloud contaminated ones. Only the NIR band of Landsat was employed as the input of C/M method because of its sensitivity to surface water.

MERIT HYDRO dataset [35,36] derived from the global hydrologically corrected DEM product (MERIT DEM), was used in this study to generate static river morphological factors, including river gradient, sinuosity, and river width. River width in the MERIT HYDRO dataset is derived through integrating flow direction and SRTM Water Body Database (SWBD) [37], which was later further improved to handle subpixel water fraction [36], enabling width estimation for rivers smaller than 90 m. It is noted that this value is a proxy of river width that does not take river channel and water level variation into consideration. In addition, we adopted water dynamic layers from Global Surface Water Dataset (GSWD) developed by European Joint Research Centre (JRC) [38] as indicators for the dynamic river morphology, which were used for qualitative analysis. Specifically, Water Occurrence and Maximum Water Extent of this dataset were used as indicators for overall river water dynamics. Water Occurrence stands for the frequency with which water was present. It shows where surface water occurred between 1984 and 2019 and provides information concerning overall water dynamics. Maximum Water Extent provides information on all the locations ever detected as water over the 36-year period.

2.2. Methodology

2.2.1. Extraction of River Morphology

In this study, river morphology indicators including river gradient, sinuosity, and width. River network is first needed before these indicators can be derived. We employed the river network data in the MERIT HYDRO dataset directly and processed the missing river network in some of the sites from MERIT DEM using the hydrology tools in ArcGIS. Based on the river centerline of the river network, river gradient (G) in this study refers to the ratio of elevation difference (H) to the horizontal distance (L). Sinuosity (S) was calculated as the ratio of the linear distance (D) to the actual river length (l). Extracting river width from MERIT DEM was following Yamazaki et al. [36]. River centerline pixels were first identified, whose corresponding river width were calculated accordingly. For each site, there are a series of river centerline pixels. Among them, the one that has the largest river width was selected, whose river width was taken as that of site. A detailed demonstration is referred to Figure A1 in the Appendix A.

2.2.2. C/M Modelling and Validation

The basic principle of C/M method is to calculate the ratio between a stable land pixel for calibration (C) and a pixel within river for measurement (M), and then apply a linear regression between C/M series and observed discharge series to give discharge information at a site [39]. In this study, time series of Landsat NIR band were used as the input of the C/M model, together with the corresponding gauge observations. As land cover is generally changing over time, to obtain a more stable C pixel, a coefficient of variation (CV) was computed and applied to each pixel, considering that CV can indicate the variability of surface reflectance [40]. Lower CV represents higher stability of land cover. Therefore, we considered those pixels that have low CV values (lower than the 5th percentile) to be stable pixels, and took their average as a surrogate of C pixel [39]. For each Region of Interest (ROI), all pixels were considered as M candidates, and a time series of C/M was calculated for each candidate. The correlation between its C/M value and the observed discharge was quantified by Pearson’s correlation coefficient (r) and validated with the corresponding p-value. Finally, the candidate pixel with the highest r value and a p-value < 0.05 (correlation is significant at the 95% confidence level) was chosen to be the optimal location of M. This ensures the C/M model built for this ROI have the best performance. In this study, each selected site near the gauge was taken as an ROI.

Once the C and M pixels are determined, a linear relationship would be established between the C/M ratio and discharge, which could then be employed to estimate river discharge from remote sensing data. Corresponding observed discharge of the validating period were used as the reference to evaluate the accuracy of C/M method in each site. Relative Root Mean Square Error (RRMSE) (Equation (1)) and Nash-Sutcliffe Efficiency coefficient (NSE) (Equation (2)) were used for accuracy assessment. In Equation (1), n stands for the number of validation data of each site, Q_rs and Q_g stand for the estimated discharge and observed discharge, respectively. $\overline{Q_g}$ is the mean value of the observed discharge. In Equation (2), the T is the total observation times, $Q_g^t$ and $Q_{rs}^t$ represent the observed value and estimated value at time t. RRMSE represents the normalized RMSE value. Lower RRMSE indicates higher accuracy. NSE is another widely used indicators for hydrological studies. Its value ranges from negative infinity to 1, with larger value representing more stable and reliable modelling results.

$$\mathrm{RRMSE}=\frac{\sqrt{\left[\frac{{{\displaystyle \sum }}_{\mathrm{t}=1}^{\mathrm{n}}{\left({\mathrm{Q}}_{\mathrm{rs}}-{ \mathrm{Q}}_{\mathrm{g}}\right)}^2}{\mathrm{n}}\right]}}{\overline{{\mathrm{Q}}_{\mathrm{g}}}},$$

(1)

$$\mathrm{NSE}=1-\frac{{{\displaystyle \sum }}_{\mathrm{t}=1}^{\mathrm{T}}{{(\mathrm{Q}}_{\mathrm{g}}^{\mathrm{t}}-{\mathrm{Q}}_{\mathrm{rs}}^{\mathrm{t}})}^2}{{{\displaystyle \sum }}_{\mathrm{t}=1}^{\mathrm{T}}{{(\mathrm{Q}}_{\mathrm{g}}^{\mathrm{t}}-\overline{{\mathrm{Q}}_{\mathrm{g}}})}^2},$$

(2)

2.2.3. Qualitative and Quantitative Analysis

Qualitative analysis is carried out by checking how C/M modelling accuracy varies with surface water dynamics, and river morphological indicators. Quantitative analysis is conducted using the Multivariate Analysis of Variance (MANOVA), a member of the General Linear Model. It is one of the most common multivariate statistical methods used for investigating whether two or more variables would significantly affect the observed variables [41]. It calculates an F value for each variable, which is taken as the ratio of its independent effect and the random effect. The independent effect is measured using Mean of Squares for factor A (MSA) which is calculated through averaging the inter-group sum of squares. The random effect is measured using Mean of Squares for Error (MSE) which is calculated through averaging the intra-group sum of squares. A significant test sig can then be derived according to the distribution of F values. A large F generally causes small sig value. If the sig is less than the significant level (usually set as 0.05), it suggests the pass of the significance test [42]. In this study, all 24 selected study sites were considered as independent samples with a series of river morphological indices (i.e., river width, sinuosity, and gradient) and accuracy indices (i.e., RRMSE and NSE), which were used as the input of MANOVA. Through comparing the F value of each river morphological indices, MANOVA is able to identify the factor that affects the discharge estimation most, with the largest F value and the least sig.

2.2.4. Demonstration of Satellite Virtual Gauge Selection

The above qualitative and quantitative analysis gives inference about how river morphology affects C/M discharge estimation. We employ another independent gauge in the MDB (Gauge 425008) to conduct a practice of applying the implications on the SVG establishment. 19 ROIs that have the same size with previously defined ROIs near the gauge are selected as candidate SVGs. We then try to find the optimal one based on the findings of the above analysis. C/M modelling and validation are carried out for all selected ROIs and the results are cross compared to confirm whether the one inferred by the river morphology is the optimal one.

3. Results

3.1. Descriptions of River Morphology

As river network in the MERIT HYDRO dataset was generated from DEM data, which sometimes embraces errors. We used river water extent visually inspected from Landsat imagery as the assistant to correct the river network. Figure 3 shows the corrected river networks in the HRB and the MDB. Accordingly, the river gradient, sinuosity, and river width of each study site were calculated and listed in

Table 2. The extraction results of gradient, sinuosity, and river width of each river section.

Gauge	Study Site	Gradient (m·km−1)	Sinuosity	River Width (m)
Qilian	1-1	15.89	1.08	64.09
	1-2	16.71	1.26	76.6
	1-3	14.31	1.10	54.99
	1-4	8.17	1.07	54.98
Zhamashike	2-1	16.2	1.06	68.52
	2-2	39.07	1.10	71.26
	2-3	7.60	1.01	55.82
	2-4	20.46	1.08	68.52
Sunan	3-1	16.87	1.00	23.4
	3-2	30.09	1.11	28.66
	3-3	22.18	1.00	63.67
	3-4	18.87	1.03	32.30
409035	4-1	0.15	1.41	62.27
	4-2	0.48	1.36	71.79
	4-3	3.57	1.11	60.92
	4-4	1.40	1.08	55.34
410130	5-1	2.12	1.28	41.15
	5-2	0.32	1.66	55.74
	5-3	9.94	1.47	53.02
	5-4	3.87	1.26	54.66
414200	6-1	0.06	1.22	78.65
	6-2	0.07	1.60	101.20
	6-3	4.48	1.11	51.78
	6-4	1.15	1.11	93.53

. The overall gradient of the HRB is larger than that of the MDB, leaving its sites generally having higher gradients that those of the MDB. For example, site 2-2 near Gauge Zhamashike has a gradient of 39.07 m/km. Meanwhile, the difference of gradient among four sites of the same gauge can be large (e.g., gradients for Gauge Sunan). As the rivers in the HRB are relatively straight, sinuosity of the MDB sites is obviously higher than that of the HRB sites. River widths of the study sites are around 60~70 m. Sunan sites have a width less than 30 m, and one of the MDB sites have a width slightly greater than 100 m.

3.2. Discharge Estimation Using C/M Method

Pearson correlation coefficient (Pearson’s r) between C/M value for each candidate M pixel and observed discharge of the 6 hydrological gauges (24 sites in total) selected in the HRB and MDB is shown in Figure 4, with red dots representing the finally selected location of the M pixel. It is observed that the M pixels tend to locate at the edge of the river channel rather than the center, probably because the pixels in the center are always wet, and their reflectance thus changes relatively small as the streamflow increases. On the contrary, pixels on the river edge would experience more intensive reflectance variation as water level changes, making them more sensitive to discharge variation. This is consistent with the findings of Tarpanelli et al. [26].

Figure 5 shows the scatter plots of the optimized C/M value and the observed discharge. The established models and their corresponding coefficient of determination (R²) are also displayed. The observed discharge and C/M series are overall positively correlated for all sites as indicated by Pearson’s r and coefficient of determination (R²), particularly for sites near Gauge Qilian and Gauge Zhamashike in the Heihe River Basin, Gauge 409035 and Gauge 414200 in the Murray Darling Basin (R² up to 0.7). However, for Gauge 410130 and Gauge Sunan, the correlation is noted to be relatively lower, which suggests considerable uncertainties may exist in discharge estimation for these sites. Moreover, four sites selected for each gauge may demonstrate quite different accuracy, with Pearson’s r differing up to 0.4 (e.g., Pearson’s r from 0.33 to 0.77 for Gauge Sunan).

Figure 6 shows the hydrograph of the estimated discharge and the observed discharge for the validation period. Estimated discharge from four selected sites and observed discharge from corresponding gauge exhibit general consistent pattern. However, when the observed discharge is high, the C/M method tends to underestimate the discharge.

Table 3. The validation results of different hydrological gauges.

Gauge	Study Site	RRMSE	NSE
Qilian	1-1	0.44	0.60
	1-2	0.54	0.39
	1-3	0.63	0.18
	1-4	0.56	0.35
Zhamashike	2-1	0.56	0.50
	2-2	0.56	0.51
	2-3	0.55	0.52
	2-4	0.55	0.53
Sunan	3-1	2.68	−10.79
	3-2	1.78	−4.21
	3-3	1.82	−4.43
	3-4	0.76	0.06
409035	4-1	0.88	−0.14
	4-2	0.60	0.47
	4-3	0.61	0.46
	4-4	0.96	−0.34
410130	5-1	1.69	−0.21
	5-2	1.59	−0.07
	5-3	2.24	−1.10
	5-4	2.00	−0.68
414200	6-1	0.65	0.40
	6-2	0.53	0.61
	6-3	1.00	−0.43
	6-4	1.02	−0.48

shows the modelling accuracy in forms of RRMSE and NSE. Overall, river discharge of Gauge Qilian, Gauge Zhamashike in the HRB and Gauge 409035 and Gauge 414200 in the MDB has been estimated relatively accurately using the C/M method, while the results for the other two gauges are less accurate. It is noted that for each gauge, its four selected sites exhibit different accuracy.

3.3. Influence of River Morphology on Discharge Estimation

3.3.1. Qualitative Analysis

Water Occurrence and Maximum Water Extent obtained from the JRC GSWD product (Figure 7) were employed for analyzing the impact of surface water dynamics on river discharge estimation. When water Occurrence within the river channel basically exceeds 90%, this means that the permanent water has been overall properly detected by the Landsat imagery. When the Maximum Water Extent is obviously larger than the channel, this means that the water dynamics have been adequately captured by historic Landsat data. Among the four river sections that have JRC data covered, study sites (3-1, 3-2, 3-3 & 3-4) near Gauge Sunan on Liyuan River have relatively lower accuracy, with RRMSE greater than 1 and NSE less than 0 (Table 3). It is observed that maximum water extent is small and water occurrence is relatively low here, with the maximum occurrence equals to 71%, and most values close to 50%. It is obvious that the river water is mostly constrained within the river channel. This suggests that surface water cannot be observed clearly by Landsat, and its dynamics is limited. This means that the water dynamics here could not be properly captured by Landsat, and thus the discharge cannot be precisely estimated, as proved in Table 3. Sites near Gauge 409035 and Gauge 410130 in the MDB are in a similar position but with a slightly higher accuracy as those near Gauge Sunan.

On the contrary, sites near Gauge 414200 have relatively higher accuracy, with some obtaining RRMSE less than 0.55 and NSE greater than 0.60 (Table 3). For this site, Water Occurrence within the river channel generally exceeds 90%, and Maximum Water Extent is obviously larger than the channel. This implies that these sites have relatively higher quality Landsat observations, which adequately reveal its extensive surface water dynamics. Therefore, we conclude that available high-quality Landsat observations and surface water dynamics are likely to be two important determinants for the accuracy of discharge estimation. Since these determinants could be inferred from Water Occurrence and Maximum Water Extent layers in the JRC GSWD, the GSWD product can serve as important support that provides prior knowledge about the suitability of SVGs.

3.3.2. Quantitative Analysis

Multivariate analysis of variance was conducted to analyze how the three river morphological indicators (i.e., gradient, sinuosity, and width) affect the accuracy of river discharge estimation. As shown in

Table 4. F statistics and Sig values of the multivariate analysis of variance, and correlation coefficient (r) between RRMSE/NSE and river morphological indicators.

	RRMSE			NSE
	F	Sig	r	F	Sig	r
gradient	0.06	0.81	0.04	2.34	0.14	−0.28
sinuosity	0.89	0.36	0.12	0.58	0.45	0.26
width	13.27	0.00	−0.56	6.32	0.02	0.54

, river width was found to be the most important factor affecting the discharge estimation, in forms of either RRMSE or NSE. It is also the only factor that passed the significance test (Sig < 0.05), which means the other two factors cannot be proved to be significantly affecting the discharge estimation accuracy. Nevertheless, it seems that the sinuosity tends to affect the estimation accuracy RRMSE (high F and low Sig value for RRMSE) and the gradient tends to affect the stability of the model NSE (high F and low Sig value for NSE). Even though river gradient and sinuosity did not pass the significance test of multivariate analysis of variance, it seems that the two indicators also affect the estimation accuracy for some of the selected study sites. Table 4 also shows the correlation coefficient (r) between RRMSE/NSE and the three river morphological indicators. It is also found that river width is the factor that has the closest relationship with discharge estimation accuracy, in forms of either RRMSE or NSE.

3.3.3. Practice on Another Independent Gauge

The above results indicate that the river width is a major factor that influences the discharge estimation when using the C/M method, and the surface water dynamics would also affect the accuracy, with high Water Occurrence and large Maximum Water Extent benefiting discharge estimation. Overall, higher estimation accuracy tends to appear at the sites with larger river width (mostly larger than 50 m). We employed another independent gauge in the MDB (Gauge 425008) to conduct a practice of applying the implications on the SVG establishment. The observed discharge of Gauge 425008 can be found in Figure 2.

19 ROIs were selected near the Gauge 425008 as candidate SVGs, with the same size as previously defined ROIs (21 × 21 pixels) (Figure 8a). These ROIs were selected based on the overall river width (>50 m) and surface water dynamics (Figure 8b,c) to ensure that the river section in the ROI is wide enough and the river water dynamics have been properly captured by historical Landsat data. River morphological indicators of these ROIs were calculated and listed in

Table 5. Gradient, sinuosity, and width of the 19 ROIs.

ROI No.	Gradient (m·km−1)	Sinuosity	River Width (m)
1	0.38	1.11	74.69
2	1.12	1.08	75.37
3	1.16	1.06	66.38
4	0.90	1.03	71.03
5	0.41	2.74	89.28
6	0.66	1.09	61.68
7	0.99	1.08	84.96
8	0.58	1.16	83.13
9	0.59	1.27	81.62
10	0.60	1.11	81.54
11	0.70	1.47	68.05
12	0.66	1.14	72.95
13	0.29	1.19	78.99
14	0.47	1.82	87.11
15	0.36	2.47	65.44
16	0.98	1.09	59.84
17	0.47	1.07	57.70
18	0.38	1.03	54.11
19	0.82	1.13	51.65

. Among them, ROI 5 has the largest river width and the highest sinuosity, a relatively higher water occurrence and large maximum water extent. According to the implications we got from the precious experiments, we infer ROI 5 to be the most suitable location as the SVG for the C/M method.

In order to further verify our inference, we calculated the Pearson’s r between the C/M ratio and the observed discharge using each pixel of the whole study site as a candidate M (Figure 9a). It is observed that pixels with high correlation (Pearson’s r > 0.7) generally fall into the previously established 19 ROIs (Figure 9b). For each ROI, we selected the pixel that has the highest Pearson’s r as the optimized M pixel, and constructed a C/M model to estimate the discharge for Gauge 425008. The accuracy of each ROI is listed in

Table 6. The accuracy of discharge estimation of 19 selected ROIs.

ROI No.	RRMSE	NSE
1	1.19	0.48
2	0.84	0.74
3	1.02	0.62
4	0.91	0.70
5	0.68	0.83
6	1.24	0.44
7	1.44	0.25
8	1.75	−0.12
9	1.88	−0.28
10	1.13	0.54
11	1.46	0.22
12	1.30	0.38
13	1.09	0.56
14	1.08	0.57
15	1.13	0.53
16	1.06	0.59
17	1.78	−0.16
18	1.04	0.60
19	0.74	0.80

, indicating that the highest accuracy occurs in ROI 5 (the estimated discharge using ROI 5 and observed discharge for Gauge 425008 is shown in Figure 9c). A chart showing the difference of the estimated discharge for different ROIs can be found in Figure A2 in the Appendix A. This suggests that ROI 5 is the best location for establishing the SVG for this river reach, which further verifies the conclusions we made from the study cases in China and Australia. Except for ROI 5, ROI 2, ROI 4, and ROI 19 showed relatively high estimation accuracy. They could be used as excellent substitutions to ROI 5, when it is continuously obscured by clouds during the modelling period. As the C/M method is an iterative optimization method whose computation cost increases drastically as the size of virtual gauge increases, quickly locate a ROI that is much smaller but generates decent accuracy would be helpful for establishing SVGs for river discharge monitoring, especially for large scale monitoring where there are massive SVGs to be established.

4. Discussion

Based on experiments conducted in a series of sites selected from China and Australia, we tried to investigate the influence of river morphology on discharge estimation using the popular C/M method. We select river gradient, river sinuosity and river width as three indicators of river morphology, note that all three indicators are static variables that do not link to river water dynamics. Even though we only selected six reaches from China and Australia considering the data availability and workload, we applied four sites to each reach and thus built 24 cases that can provide fair support to investigate the scientific questions we aim to answer. We identified that river width is the most import factor that controls the reliability of discharge estimation using the C/M method. It is also noted that the other river morphology indicators, such as river gradient and sinuosity, may also affect the accuracy of the C/M method. This finding is overall consistent with existing studies (e.g., [12,31,32]). Even though they focused on different discharge estimation methods, there is a consensus that river morphology is critical for discharge estimation, and therefore should be taken into consideration when selecting SVG reaches.

River water dynamics, as constrained by the river morphology, were also found to be an important factor that may affect the C/M accuracy. As the river water dynamics can be easily obtained from remote sensing imagery and even readily available in several surface water products (e.g., GSWD), they can serve as useful resource for preselecting potential suitable SVG locations. The finding of this study was further applied an independent river reach. It was proved that it can help quickly locate a smaller ROI that generates the best C/M model for discharge estimation. Considering that the C/M method is a iterative algorithm whose computation cost increases drastically with the size of ROI, this study provides useful practice guidance for establishing SVGs for river discharge monitoring, especially for large scale monitoring where many SVGs are needed.

There are several issues that may limit our conclusions or need further attentions. First, extracting river morphology from remote sensing images and DEM may introduce uncertainties, especially in arid and semi-arid areas, where the river is either meandering or relatively narrow. Second, gradient and sinuosity failed to pass the significance test of multivariate analysis of variance, probably due to the limited number of cases. Third, water and land color change will affect the discharge estimation, as the C/M method utilizes the ratio of NIR reflectance of a C pixel and a M pixel [43]. In order to weaken the impact of land cover change (in particular vegetation phenology), we adopted the CV value of each pixel along the time series, and considered the 5th percentile of low CV pixels as land pixels. These pixels were averaged to represent the C pixel. Through this process, it is believed that those land pixels whose NIR reflectance changes drastically will not be used for the C/M modelling. M pixel was determined through an iterative optimization process. Watercolor, in particular the turbidity, would certainly affect the NIR reflectance, with higher turbidity usually causing higher NIR reflectance. Therefore, the M pixel embraces different levels of uncertainty according to the turbidity level in the SVG. Higher turbidity makes the M value greater than it should be when water is clean, and thus may lead to model overestimation. One thing that further complicates the situation is the uncertain and unstable relationship between turbidity and streamflow. Therefore, in addition to river morphology, watercolor may be another important factor that affects the spectral-based discharge estimation method. The SVG site selection may need to consider the watercolor variation as well, in particular for those rivers whose turbidity varies fiercely across space and time.

We only adopted the C/M method in this study. There are several different discharge estimation methods that are based on quite different principles. Whether the river morphology also affect the accuracy of these methods is still not clear. In the future, the input of experimental river sections, observed discharge data, and remote sensing images can be increased, and multi-source remote sensing images may also be integrated into different discharge estimation models, to comprehensively investigate how river morphology affects different approaches of discharge estimation.

5. Conclusions

When establishing SVGs and implementing discharge monitoring in the basin, finding the optimized location that is able to provide accurate river hydrological parameters for estimating discharge is critical for improving the efficiency of models. Combined with surface water dynamics derived from Landsat imagery and static river morphological indicators derived from DEM data, this paper qualitatively and quantitatively investigated the influence of river morphology on the C/M discharge estimation model, especially for medium to small rivers (width < 100 m). Through experimenting at 24 sites in semi-arid regions of China and Australia, guidance for selecting appropriate river reaches to establish SVGs in another river section (near Gauge 425008) were provided. Conclusions of this study are as follows.

(1)

The observed discharge and the C/M series in all the selected river sections showed medium-to-high positive correlation. For each gauge, the estimated results of its four study sites exhibit a moderate but different accuracy. This suggests that the C/M method is overall reliable for all the selected study sites.

(2)

Estimation results are generally better in study sites with higher water occurrence and larger maximum water extent. This implies river sections with high-quality Landsat observations and extensive surface water dynamics are important prerequisites for reliable discharge estimation. Therefore, some surface water products (e.g., GSWD) could be useful resources for assisting SVG location.

(3)

It is found that river width has the most significant influence on the discharge estimation accuracy and model stability (sig < 0.05). Even though river gradient and sinuosity did not pass the significance test of multivariate analysis of variance, it seems that the two indicators also affect the estimation accuracy for some of the selected river sections. This suggests that river morphology is truly an important factor that influencing the C/M river discharge estimation method.

(4)

Through testing at another independent river reach, our findings about how river morphology affects discharge estimation and how can they be applied for establishing SVGs are verified, which further confirmed the significance of this study. It provides useful practice guidance for establishing SVGs for river discharge monitoring, especially for large scale monitoring where many SVGs are needed.