remote sensing A New Method for Long-Term River Discharge Estimation of Small- and Medium-Scale Rivers by Using Multisource Remote Sensing and RSHS: Application and Validation

: River discharge is an important hydrological parameter of river water resources. Especially in small- and medium-scale rivers, data deﬁciency is the biggest problem for studies of river discharge. In recent years, remote sensing has become a rapid and convenient method to estimate river discharge. However, remote sensing images still have some difﬁculty generating continuous long-term river discharge. To address this problem, we developed a new method coupling the remote sensing hydrology station method (RSHS) with statistical regression downscaling, using data from optical satellites (Landsat-8, Sentinel-2), radar satellites (Sentinel-1), and un-manned aerial vehicles (UAVs). We applied this method to monitor monthly river discharge for small- and medium-scale rivers from 2016 to 2020 on Yunnan-Guizhou Plateau and evaluated the accuracy of the results. The results show that (1) by applying the newly constructed method, the water body continuity index obtained by Landsat-8 increased by 7% and the average river length percentage in the channel reached 90.7%, a 40% increase; (2) there were only 10 river ﬂow data points, on average, in the 5-year period obtained before this method was applied; after this method was applied, more than 50 river ﬂow data points could be obtained, on average, extending the quantity of data ﬁvefold; in addition, improper extreme values could also be avoided; (3) with better continuity of water body distribution, the images provided steadier river widths. The relative error of daily ﬂow estimation from Landsat-8 images was reduced by 60% and the mean percentage error was reduced by one-fourth. The relative error of the multisource remote sensing composited ﬂow was reduced by 37% with a reduction in the mean percentage error of over a half; (4) in addition, we found that when the threshold difference between water bodies and land in remote sensing images is more than 0.2, the impact of water body recognition error on ﬂow accuracy can be ignored. This method helps to overcome the absence of remote sensing methods for the long-term estimation of ﬂow series in small- and medium-scale rivers, improves the accuracy of remote sensing methods for calculating ﬂow, and provides ideas for regional water resource management and utilization.


Introduction
River discharge is an important index of river water resources. Long-term and continuous river discharge data are of great significance for watershed water resource management and utilization, and the economic development of ecosystem services [1,2]. At present, there are three principal ways of monitoring long-term river discharge: hydrologic station measurement, hydrologic models, and estimation using remote sensing data. Traditional water bodies. Therefore, the calculated river discharge results had difficulty meeting the requirements of long-term continuity. In view of the lack of definition of remote sensing images, some researchers have attempted to solve the problem. For example, for images completely blocked by clouds and fog, many scholars have developed cloud removal algorithms [19][20][21], although, due to the absence of the infrared band on Sentinel-2 and other satellites with 10 m resolution, the cloud removal algorithm is difficult to apply in this series [19,22]. Currently, the cloud removal algorithm can only be used for images covered by clouds on a large scale, but for middle-and low-latitude areas such as the Yunnan-Guizhou Plateau, the low-definition improvement in occlusion caused by thin clouds and mist on a small scale is limited.
Further addressing images with low use efficiency, we can turn to other data with lower cloud content to improve its clarity and increase the use efficiency of the images. In the fields of hydrology and meteorology, statistical linear regression, a commonly used method, is used to build a relationship between the two kinds of precipitation products to downscale and reduce error, which helps increase the data use efficiency [23,24]. The relationship between coarse resolution precipitation data, with a long time span, and fine resolution precipitation data, with a short time span, can be established to apply coarse resolution precipitation data on a smaller scale and retain its characteristics over a long time span. This method has also been widely used to improve the utilization efficiency of remote sensing images and there have been many studies to build mapping relationships between satellite images with different resolutions. By using relatively fine-resolution images such as the Landsat series, coarse-resolution images such as those from the Moderateresolution Imaging Spectroradiometer (MODIS) are downscaled [25][26][27]. This method is simple in principle and can improve the utilization efficiency of remote sensing data that originally had a long time span but had problems with clarity. Moreover, it also helps process discontinuous fragments into continuous images while retaining the long time span, thereby helping the long-term flow estimation of small-and medium-scale rivers.
In addition to optical remote sensing, radar remote sensing can also become an important data source for calculating river flow. The Sentinel-1 synthetic aperture radar (SAR) satellite, launched in 2016, has a high spatial resolution of 10 m. At the same time, as an active remote sensing instrument, it is not affected by clouds, tree canopies or time of day. It offers new options for solving the problem of discontinuous monitoring of small-and medium-scale rivers and removing the interference of clouds and fog [28], and it can be combined with traditional optical satellite data to improve the accuracy of remote sensing image flow estimation.
In summary, remote sensing, as an important means of estimating river discharge, is not restricted by time and space in the same way as the traditional station model. Combined with NDWI and other water body indices, water body information can be accurately obtained. However, because of the limitation of the imaging principle, problems such as low image resolution and low use efficiency of images are encountered in many low latitude regions. This paper tries to couple statistical linear regression and the Remote Sensing Hydrology Station method, integrate many kinds of remote sensing data, and build a small-and medium-scale river flow calculation method. In view of the problem of low resolution images, it helps improve efficiency and generate continuous long-term sequences. Additionally, the accuracy of the calculation results is verified.

Study Area
This study selected the Xingjiang River basin in the Yunnan-Guizhou Plateau region in southwest China as a typical research area ( Figure 1). The Xingjiang River is a tributary of the Yangtze River-Wujiang River system located in the upper reaches of the Hongfeng Lake Reservoir. It is an important supply source of the Hongfeng Lake Reservoir, which is an important drinking water resource for Guiyang [29]. Rainfall is abundant in the basin, up to 1190 mm annually. There is a saying in this area that There Were Never Three Sunny Days in a Row [30]. Because the basin is located deep inland and in the middle of the Yunnan-Guizhou Plateau, there is only one field survey station, Huangmaocun National Hydrology Station, within the basin, classifying it as a small basin with insufficient data. Abundant clouds in this area, especially perennial thin clouds and fog, result in low clarity of satellite remote sensing images, which greatly affects the utilization efficiency of remote sensing data [31]. It is an important water source area in the upper reaches of the Yangtze River Economic Belt, adjacent to the Pearl River System. This basin is thus an important economic and cultural development area where two river systems meet, and river water is a vital ecological resource in this area.
Remote Sens. 2022, 13, x FOR PEER REVIEW 4 of 20 Lake Reservoir. It is an important supply source of the Hongfeng Lake Reservoir, which is an important drinking water resource for Guiyang [29]. Rainfall is abundant in the basin, up to 1190 mm annually. There is a saying in this area that There Were Never Three Sunny Days in a Row [30]. Because the basin is located deep inland and in the middle of the Yunnan-Guizhou Plateau, there is only one field survey station, Huangmaocun National Hydrology Station, within the basin, classifying it as a small basin with insufficient data. Abundant clouds in this area, especially perennial thin clouds and fog, result in low clarity of satellite remote sensing images, which greatly affects the utilization efficiency of remote sensing data [31]. It is an important water source area in the upper reaches of the Yangtze River Economic Belt, adjacent to the Pearl River System. This basin is thus an important economic and cultural development area where two river systems meet, and river water is a vital ecological resource in this area.  We selected segments 0501, 0502, 0401, 0402, 0404, 0405 and 0406 through field investigation and comparison, as shown in Figure 1, comprising of one natural straight channel, two natural curved channels, two straight channels affected by artificial influence, and two curved channels affected by artificial influence. Segment 0406 is close to Huangmaocun National Hydrologic Station, which was used as the in situ data site for comparative analysis. These segments cover the main distribution region of the Xingjiang River system and contain the key positions where tributaries are added. These segments can therefore effectively monitor the discharge changes of the Xingjiang River and help verify the accuracy of the flow values calculated by the application of this method.

Data
This research adopts as data sources three kinds of satellite images and unmanned aerial vehicle (UAV) aerial images, along with the observed hydrological data. In terms of SAR data, we chose Sentinel-1A as the supplement images for their high spatial resolution and free access. By using a DJI Mavic 2 UAV with an HD camera and GPS positioning, the river and the surrounding environmental conditions in each study area were scanned with high spatial resolution (5 cm × 5 cm). With the help of the professional UAV image processing software Pix4D, we were able to obtain high-precision images after image synthesis, point cloud encryption, and spatial information calibration to generate centimeter-level digital surface models (DSMs) and digital orthographic images (DOMs). We also gathered information on the location and elevation of the segment of the digital river model. At the same time, open and free satellite images, including Landsat-8 (LS8), Sentinel-1 (ST1), and Sentinel-2 (ST2) series images, were used as data sources to identify the distribution of river channels and water bodies. With the help of Google Earth Engine (GEE, https://earthengine.google.com/ (accessed on 21 February 2022)), image projection, image stitching, cutting and radiation correction calibration were rapidly finished in the cloud platform, and the green band and near infrared band were computed for the normalized index of the water body figure (NDWI) online. For the SAR data, we processed the back-scatter value into the range from −1 to 1, and chose the typical area of water body and land, and used these to obtain river reach data. In this case, we could identify water body range from SAR images as a supplement to optical images. The flow velocity and water depth were measured on site in each section. The flow velocity was measured by an LS300-A rotor current meter. The measured flow was calculated by the three-point method and processed as an average flow velocity in each section. The water depth was measured by a Deeper Sonar detector. The average water depth of multiple points was combined with DSM river elevation information obtained by UAV to simulate the channel water shape distribution and cross-section area, which was multiplied by the measured flow rate to obtain the channel flow rate. For the in situ measured data, we selected the monthly discharge values from 2016 that coincided with the time range of the ST1, ST2 and LS8 satellites. The specific in situ discharge values were extracted from the daily discharge table of the Yangtze River System (Wujiang River Basin) in the Hydrological Yearbook of 2016. Table 1 summarizes these datasets. We used the Remote Sensing Hydrology Station method [6,17,[32][33][34][35] to calculate river flow in the study area, which included three main steps: construction of hydraulic geometry, relationship between flow and parameters, and key parameter acquisition and volume estimation.
Step 1: Before beginning the work, a steady and typical segment should be ensured. We obey the following rules to confirm the target segment: 1 steady river banks, 2 straight and long-enough channel, 3 few constructions and canopies, 4 single channel. Construction of hydraulic geometry is carried out by Mavic 2, a consumer-grade UAV produced by DJI Company (https://www.dji.com/cn (accessed on 21 February 2022)). Currently, there are many studies that use UAV ( [36][37][38]) to carry out low-altitude flights over river sections that need to be calculated, and combine high-precision ground images and DSMs to obtain elevation information. Meanwhile, on the ground, real-time kinematics (RTK) are applied to carry out field measurements and combined with UAV DSM to determine the specific values of water depth, velocity, slope, and roughness. By integrating these data, the river section shape can be fitted.
Step 2: The Manning formula, a classical empirical formula in hydraulics, provides a simple and highly accurate method for estimating flow [39][40][41]. The basic formula for calculating river flow is Equation (1). The main influencing parameters are roughness, hydraulic gradient, area of water passage and wet-cycle length. According to the natural conditions of the river channel, it was categorized into one of three types: triangle, trapezoidal or rectangle [6,39,42,43]. The roughness and hydraulic slope were confirmed from our in situ measurement and UAV low-altitude images, which would not change in short time periods. Based on the measured cross-section data, the area of water passage and wet-cycle length were ensured with the typical relationship of the channel shape. In this case, the changing curve of the single relationship between the flow and river width was calibrated.
In the equation, k is the conversion coefficient, where 1 is taken, n is roughness, A is area of water passing section, P is wet-cycle length, and S is the hydraulic slope, which is replaced by the channel slope.
Step 3: After the construction of the digital river channel, the river width became the only parameter needed for calculating the discharge of the long series. The historical remote sensing images of the station were mainly obtained through GEE and calculated into NDWI values, which were used to identify the water distribution in the river channel and select water from the images. After the land training area and water training area were selected, the water body threshold was obtained in batches then used to calculate the width of water pixels in the channel in the training area. For those pixels within these thresholds, we took a simple linear equation to obtain the ratio of water and land. As a result, these pixels were divided into minor water pixels and land pixels. For average condition of the segment, we calculated the total water area of the segment and the length of its middle line to get the average width of water range. Finally, changes of water body width in the river and non-water riverbank were monitored. Table 2 summarized the parameters mentioned in RSHS method and their sources.

Parameter Symbol Source
River width -NDWI calculated from satellite images Area of water passage A River width and depth measured in situ with the mathematical relationship based on the shape of river channel Image downscaling refers to processing remote sensing images with large scale and low spatial resolution into images with small scale and high spatial resolution, which can reflect the finer differences of spatial information in images. The principle of statistical downscaling that was applied in this method is to build a correlation between different resolution datasets. To deal with the LS8 image with 30 m resolution, the images were processed by statistical linear regression into the same 10 m resolution as ST1 and ST2. The specific method was as follows: First, LS8 images were resampled from 30 m resolution to 10 m resolution, and the resampled original image (LS8 resampling image) and the target image (ST2) were aligned in the image grid range. Then, the linear regression method was applied to both the original and target images from the same period to obtain the correlation coefficient. Finally, according to the difference between the target image and the calculated value, the residual image was obtained, which was taken into the original image for linear calculation to obtain the regressed image. The specific equations are as follows: In the equations, Obj is the target image (ST2, ST1 image). Ori is the original image (LS8 resampled image). Coe1 is the correlation coefficient between LS8 and ST1. Coe2 is the correlation coefficient between LS8 and ST2. Res is the residual image. Reg is the regression image. i and j are the row and sequence numbers of grid points in the image, respectively.

Water Identification Assessment: Water Continuity Index & River Length Ratio
Water body recognition by remote sensing has been a very important means of water body monitoring in recent years. Generally, the accuracy of water body recognition by remote sensing is evaluated through indicators such as river width, water area and river continuity. We generated high-precision images by UAV at the test site to study the extent of the river and extract the NDWI threshold in the training area of the water body according to the method above. The pixel extent of the water was calculated to obtain the average width of the river. Then, the average width was used to divide the water area to obtain the length of the water in the channel and to calculate the relative river length within the whole river, which is named the river length percentage (L i ) in this paper. In addition, we define an index C to evaluate river continuity based on the landscape fragmentation index [28].
Remote Sens. 2022, 14, 1798 8 of 19 In the equations, i is the number of the corresponding section, N i is the number of continuous water patches at section i, A i is the total number of water grid points at section i, RW i is the river width, and LH i is the overall length of the river.

Runoff Accuracy Assessment: Error Index
To evaluate the rationality and reliability of this method for runoff estimation, root mean square error (RMSE) and mean percentage error (MPE) were selected as precision evaluation methods. RMSE was used to analyze the overall reliability. MPE was used to evaluate the difference between the estimated flow and the measured flow of a single section. Equations for these methods are as follows: In the equations, Q e is the estimated flow, Q m is the measured flow, and n is the calculation number.

Result of Linear Regression Downscaling Method
For each of the seven segments, we applied identification of water bodies and riverbanks on the LS8 original images and LS8 regression images. The results showed that there were significant differences between the two images: in the original images, water bodies were discontinuous, the distribution of river width was uneven, and there were fewer pixels in the channels. In regressed images, the water bodies were more continuous, river width was evenly distributed, and the number of water pixels increased significantly, as shown in Figure 2.
Among the seven segments, the water bodies identified by LS8 original images were discontinuous, with an average water body continuity index of 0.918, while the index of regression images reached 0.983, an increase of 7.1%. Affected by the spectral characteristics of the water body, NDWI images from optical satellites, such as LS8, recognize regions with deep water depth, slow flow, and smooth water surface as water bodies, as shown in Figure 2(b1,c1,d1,g1). In these areas, the river channels are wide and there are bends in the river. Deep pools are can easily be found. In addition, the reaches with artificial dams in the river and high obstruction of water bodies have higher NDWI values, which are easy to identify [45]. Many of the areas identified as land in the river channels have shallow water depths and fast flow velocities. The NDWI value of the region with an unstable water surface is lower, so the water pixels of these reaches show poor water continuity and uneven river width distribution.
In Figure 2(a1,e1,f1), rivers are narrow or straight, and the water depth is not deep. Therefore, the water bodies identified from LS8 original images are mainly distributed in the center of the river, where the water pixels show uneven river width distribution and low continuity. In addition, the resolution of the LS8 original image is coarser than that of ST2, at approximately 30 m. The pixels close to the riverbank are usually in a mixed state of land and water, and it is difficult to separate them under the resolution of the original LS8 image. As seen from Figure 2(h1,h2), the average river length percentage in the regressed images increased from 51.2% to 90.7% after the method was applied, an increase of 40%. The LS8 image resolution reached 10 m after method application, with the help of ST2 images in the same period to determine the location of the water in the river. The NDWI value of areas with shallow depths or unstable surfaces also reached the level of other water bodies. The scope of identified water expanded, with water pixels distributed more evenly in the river. In addition, the pixels that originally mixed land and water were divided into water pixels with high NDWI values and land pixels with low NDWI values. The remaining water bodies were connected to the originally recognized water pixels, which improved continuity and solved the problem of water bodies with uneven distributions of river width being less recognizable due to the resolution effects.
LS8 image. As seen from Figure 2(h1,h2), the average river length percentage in the regressed images increased from 51.2% to 90.7% after the method was applied, an increase of 40%. The LS8 image resolution reached 10 m after method application, with the help of ST2 images in the same period to determine the location of the water in the river. The NDWI value of areas with shallow depths or unstable surfaces also reached the level of other water bodies. The scope of identified water expanded, with water pixels distributed more evenly in the river. In addition, the pixels that originally mixed land and water were divided into water pixels with high NDWI values and land pixels with low NDWI values. The remaining water bodies were connected to the originally recognized water pixels, which improved continuity and solved the problem of water bodies with uneven distributions of river width being less recognizable due to the resolution effects.  Figures (a-g) represent the water body distribution in the map before and after the method was applied to segments 0501, 0502, 0401, 0402, 0404, 0405 and 0406. In each group, 1 represents the original image, 2 represents the regressed image. Group h represents the statistical variation of seven segments in two parameters: continuity (h1) and river length percentage (h2). Bars in light blue represent original images, bars in dark blue represent regressed images.

Effect of Extension to Runoff Series from Multisource Data
We constructed the river flow sequence of ST2 optical images, ST1 radar images, LS8 original optical images and LS8 regressed images and integrated the flow data into a com-posited flow sequence according to the in situ value series of typical years. The monthly flow data series of seven sections from 2016 to 2020 were obtained, as shown in Figure 3. As indicated in Figure 3a-g, only approximately 10 flow values were estimated in the sequences taken from the original optical images (LS8), while there were more than 50 comprehensive sequence values in the composited sequence, five times the original data. Compared with the results of various data sources, the composited sequence greatly extended the timespan of optical images and ensured the temporal continuity of flow.  As seen from the figure, ST2 optical images and LS8 original optical images, the traditional data sources used by flow calculation, are affected by a large amount of cloud and fog in the research area. For 60 months over five years, the available data amount from seven sites was very small, with fewer than 30 data points from ST2 images and fewer than 20 data points from LS8 original images. In addition, due to thin cloud occlusion, the NDWI values of most pixels were generally higher when synthesizing optical images. This led to the overestimation of water body recognition when the recognition algorithm was learning to recognize water body or land pixels and mistakenly identified many areas covered by clouds as water bodies or land. As a result, the calculated flow reached improper extreme values within 0.1 m 3 /s or more than 100 m 3 /s.
As an active remote sensing data source, ST1 radar images are not affected by weather conditions, so without special interference, ST1 images can cover the timespan of the whole year and can be used as a complete flow calculation data source. However, radar satellites are different from optical satellites, as their imaging principle is to send electromagnetic pulses and receive echoes in the sensor, calculating underlying surface roughness based on the strength of the returned signal (backscatter). Backscatter from water features is similar to some other surface features, such as artificial ground, which reduces identification capability for water compared to optical images, especially in NDWI images. There is also often noise in water body recognition, as shown in the figure, where the non-flood season flow shows outliers greater than 20 m 3 /s. The flow of LS8 regressed images after downscaling is within the flow range of typical years, which can be used as a reference value for correction of flow calculated from radar images in the non-flood season. The composited sequence retains the radar data long time span and excludes the improper extreme values in radar data with the aid of ST2 images with LS8 regressed images, as shown in Figure 3 in the sequence (a), (b), (c), (d). Reasonable optical image data are used as a supplement to the composited sequence to help reduce the emergence of improper extreme values in the flow data. If there are few improper extreme values in each data source, the composited sequence will filter and adjust the reasonable estimated flow results of each data source according to the annual distribution of typical years and synthesize the flow sequence that meets the annual distribution law of true values.

Accuracy Comparison between Runoff Results from Original and Regressed Images
After statistical regression downscaling of LS8 original images in a typical year (2016), it can be seen that the daily calculated flow results were closer to the measured true value, and the RMSE and MPE were greatly reduced, with a decrease in a range from 3 m 3 /s to 8 m 3 /s. It is shown in Figure 2 that the continuity of water bodies improved after the method was applied. As seen from Figure 4 a,b, the daily flow accuracy of LS8 after regression downscaling was significantly improved. RMSE before regression was 13.5 m 3 /s, and after regression, it decreased to 5.4 m 3 /s, reducing by nearly 60%. This indicates that the discontinuity of the original images provided unreliable average river widths, which led to unstable discharges with an error of 13.5 m 3 /s on average, and this number decreased greatly with steadier average river widths provided by water bodies with better continuity. The MPE was reduced from approximately 216% to 44%, which shows that the average difference between the calculated discharge and observed discharge decreased from 10.9 m 3 /s to 4.22 m 3 /s. In the flood season, the RMSE of the two flow points decreased from 8.3 m 3 /s to 6.2 m 3 /s, with a reduction of 25%. The MPE decreased from 59% to 44%, indicating that the average difference between calculated discharge and observed discharge decreased from 6.73 m 3 /s to 3.27 m 3 /s. According to the distribution of flow points, flows estimated before regression come with obvious random errors, which are distributed on both sides of the 1:1 line and in general relatively far away from the 1:1 line. After downscaling, the estimated flow was distributed near the 1:1 line. The whole fitting line slope was close to 1, with a smaller random error, and the fitting line R 2 was 0.77. This shows that the fitting result of the flow point is reliable.   The accuracy of the composited monthly flow series also improved after regression. As shown in Figure 4c,d, the overall RMSE of monthly flow reached 11.7 m 3 /s before regression but decreased to 7.4 m 3 /s after regression, a decrease of nearly 37%. At the same time, MPE was reduced from 213% to 107%. This means the average difference between the calculated discharge and observed discharge decreased from 10.0 m 3 /s to 5.99 m 3 /s which is a significant effect. In the flood season, the RMSE of the three flows decreased from 9.43 m 3 /s to 5.71 m 3 /s. The MPE decreased from 46% to 29%, as the average difference between the calculated discharge and observed discharge in flood season decreased from 8.93 m 3 /s to 4.27 m 3 /s. In addition, from the distribution of flow points in the figure, the accuracy was significantly improved: the systematic error of the original monthly flow is obvious, and more flow points are distributed on the right of the 1:1 line, indicating that the estimated flow is too large. The fitting line R 2 is only 0.008, indicating that the data are relatively discrete. After regression, the systematic error was significantly reduced so that the flow points were closer to the 1:1 line, with the slope of the fitted line closer to 1, and R 2 reaching approximately 0.4, which indicates higher concentration and more reliable results.
It can be found from the figure that the accuracy of flow in the flood season is generally higher than that in the non-flood season, and the systematic error is smaller. In the nonflood season, the estimated flow is more likely to be larger, which is discussed in Section 4.3.

What Affects the Accuracy of This Method?
As shown in Figure 3, under the influence of thin cloud and noise errors in the study area, the water clarity of the image is low and the flow estimation is far overestimated. When the water clarity and natural conditions in the images change, the water area distribution and flow accuracy estimated by the method show significant changes. In fact, water clarity indicates the difference between the reflection characteristics of water and land, which can be observed most directly from their NDWI value. As we have chosen different water and land area in the images, the difference between these thresholds varies. This study defines water clarity alone as the degree of distinction between the NDWI threshold of water pixels and land pixels in the remote sensing images, which is named P i to be the judgment standard.
In the equation, P i is the water clarity of each segment, V water is the NDWI threshold value of water, and V land is the NDWI threshold value of land.
We analyzed section 0406 and set up seven situations in this section. In each situation, the NDWI thresholds of water bodies and land surfaces were different due to different selections of water bodies and land surface areas. In Figure 5a, water clarity values are arranged by gradients from 0 to 0.25. As the water clarity value increases, the flow error is reduced. Its change tendency is not always the same, especially in situations where the water clarity value is less than 0.05 and the flow relative errors are 100% or more. In these situations, a slight increase in the water clarity value leads to a significant decrease in the flow rate error, which is due to the similarity of water and land NDWI values, resulting in a large water body identification error. As the images in the water become clearer and the water clarity value increases to over 0.05, the surfaces of the water and land differ more and the flow error starts to fall. When the water clarity value is over 0.2, the flow errors can reach 6% and the rate of decrease slows. This indicates that when water clarity exceeds 0.2, the difference between the water body identification error and the flow estimation accuracy is negligible, and the influence of the flow rate on estimating the effect is improved. the water clarity value increases to over 0.05, the surfaces of the water and land differ more and the flow error starts to fall. When the water clarity value is over 0.2, the flow errors can reach 6% and the rate of decrease slows. This indicates that when water clarity exceeds 0.2, the difference between the water body identification error and the flow estimation accuracy is negligible, and the influence of the flow rate on estimating the effect is improved.  Figure 5. (a) Change diagram of flow relative errors and water clarity values, in which the red line represents the water clarity value and the blue line represents the relative error of the flow rate; (b) Statistical graph of the water clarity values of each section, in which the light blue column represents the original water clarity value of each section, the dark blue column represents the regressed water clarity value of each section, the orange column represents the original value of water in each section, the red column represents the regressed value of water in each section, the light green column represents the original value of land in each section, the dark green column represents the regressed value of land in each section, and the red dotted line represents the dividing line where the water clarity value is 0.2.
As seen from Figure 5b, in the original sections where the regression method was not applied, most of the water clarity values are below 0.2, with only two sections exceeding 0.2. At present, the downscaling method we adopt can significantly raise the water clarity values up to above 0.2 in the corresponding sections, but there are still some sections whose water clarity values show little effect from the regression method, due to environmental conditions or cloud cover. According to the relationships in Figure 5, the applied method in this paper helps raise water clarity values up to above 0.2, and the flow errors in the corresponding sections can thereby show a substantial decline. Further application in the future and improvement of this method should focus more on adjusting and controlling water clarity values in the range above 0.2, verifying the effect of water clarity values on the reduction of flow errors.

What Are the Main Influences on the Long-Term Discharge Generation of Small-and Medium-Scale Rivers?
The main factors affecting the production of long-term series flow data for small-and medium-scale rivers in this region are meteorological conditions such as clouds and fog, geographical characteristics such as river width distribution, noise errors of optical and radar satellite sensors, and the availability of remote sensing image data.
When remote sensing images are taken to directly estimate small-and medium-scale river flows at middle and low latitudes, the most influential factor is that these areas are seriously affected by clouds, shadows, and fog. Even relatively thin clouds and fog, due to their various types and different spectral attributes, cause different changes in remote sensing image values, such as NDWI [18,46]. The value of cloud shadow is close to that of ground objects such as water bodies, which leads to the error of misidentifying water bodies even in images without direct occlusion by clouds. [19]. Most ST2 and LS8 optical images in the study area are shielded by clouds, making it difficult to see the distribution of river channels and water bodies. This also leads to the fact that even though LS8 and ST2 satellites have a short revisit cycle, with averages of 25 and 47 images per year, respectively, only four and nine of these image sets could clearly identify river channels each year, which made it difficult to meet the requirements for long sequence inversion.
Small-and medium-scale rivers are more difficult to accurately identify due to their narrow channels and small water surface area, which were not easy to obtain from images due to resolution limitation and interference from other ground objects. As can be seen from Figure 2, many pixels near the river banks were recognized as land before our method was applied. This is because that most of these pixels were mixed pixels with both water and land. Under Landsat-8's spatial resolution of 30 m, these pixels could cover half of the river area and caused large errors. Under high resolution of 10 m with our method applied, more pixels were divided and many of them were classified as water. Meanwhile, there were still some pixels with mixed water and land, which is inevitable. Compared with the original results, more water area has already been correctly identified.
Moreover, even with the improvement in spatial resolution, a large number of artificial structures, such as paddy fields and roads around river channels, are also mistakenly identified as water areas, which causes more interference with river channel recognition [44,45,47], resulting in extreme outliers in water body recognition in the images and discontinuous estimation of river flow. In the period not covered by clouds and fog, the river in the image is only shown as one pixel wide in the LS8 images. In areas with a high degree of vegetation coverage on the bank or where a bridge spans the river, the edge of the water body cannot be identified accurately, so the direction and width of the flow cannot be extracted from the remote sensing images. In the ST2 image, the water body performance is better, and the water area identification and flow estimation results of the river center are improved due to the higher resolution, but there are also problems due to vegetation and artificial structures.
In addition, noise errors caused by the imaging principles of different sensors also affect the production of long series flow data in small-and medium-scale rivers. Taking ST2 and ST1 as examples, they represent passive optical remote sensing sensors and active radar remote sensing sensors, respectively. As mentioned above, the noise of optical images mainly comes from clouds and fog, while the noise of radar SAR images mainly comes from the shadows of mountains and buildings, flat road surfaces and ships in water [28,34,44]. A smooth road surface around the river or the shadow of a mountain may make the identification of the water area too large or too small, which does not meet the conditions of river calculation and has a certain impact on the calculation of long series flow.
It should be noted that at present, most spatial resolutions of free remote sensing data are too coarse for small-and medium-scale rivers, and satellite sensors such as MODIS are completely unable to be used at resolutions below 250 m. The finest Sentinel series only started in 2016 and access to the more sophisticated GF-2 commercial satellite is expensive and it has not seen widespread use. Therefore, LS8 is still mainly used for long sequence flow production.

What Are the Differences between Discharges in Flood and Nonflood Seasons?
In part 3.3, we carry out the accuracy verification of the daily and monthly flow sequences of section 0406 in 2016, as shown in Figure 4. We find that the flow estimation accuracy was higher in the flood season than in the non-flood season. This is because the actual channel shape distribution was more complex. As we adopted the method of hydraulic geometric fitting of the river, there were some errors in the fitted shape of the channel under the water surface. The underwater part may be controlled by different hydraulic geometry relations corresponding to the wet season and dry season. Here, we calculated the river width error and flow error according to the true and estimated flow values of section 0406 in 2016 and their corresponding river widths, which were grouped according to the flood season and non-flood season and ranked according to the river width error from small to large. Figure 6 shows that the flow rate error growth rate was far lower than that of the river width error in the flood season, while the error flow growth rate during most of the non-flood season was higher than that of river width, which indicates that error caused by river width recognition had a larger impact on flow accuracy in the non-flood season than in the flood season. Our algorithm is more sensitive to large flow calculations; relatively speaking, the flow calculation error of low discharge may be larger. This is because when the flow is low in the non-flood season, the water was mainly distributed in the middle of the channel, and the shape of the center of the channel was deep and narrow with a V-shape rather than wide and shallow with a U-shape. With the slight change in the river width in this period, hydraulic aspects of the change in the flow were relatively more significant. In the flood season, the change in discharge with river width was not significant. The influence of complex river shape on the relationship between river width and discharge exists not only with in a single section but also with a certain deviation between upper and lower sections of the same reach. The hydraulic geometric relationship of stations may also change with changes in time, climate and topographic conditions. If the river shape and hydraulic geometry are fitted according to the average width of the flood and non-flood seasons, and different hydraulic formulas are used to calculate them, the complex morphological changes of the river channel can be taken into account, and the errors associated with the flood and non-flood seasons can be reduced.

Conclusions
To solve the problems associated with long-term small-and medium-scale river flow monitoring, we developed a new method based on statistical linear regression and Remote Sensing Hydrologic Stations to achieve long sequence flow from high-precision remote sensing images. With remote sensing images downscaled from 30 m to 10 m resolution and multisource data integration of optical satellites, radar satellites and UAVs, river flow can be calculated with high accuracy and composited according to the rationality of flow data. In this case, different information can be obtained from multi-source remote sensing images, which indeed raises the utilization efficiency of these data sources. The method was tested on the Yunnan-Guizhou Plateau in southwest China and the accuracy of the results was verified. The results showed that (1) after this method was applied, the spatial resolution of Landsat-8 was significantly improved by linear downscaling and the water continuity index was increased by 7.1%. The average river length percentage in the

Conclusions
To solve the problems associated with long-term small-and medium-scale river flow monitoring, we developed a new method based on statistical linear regression and Remote Sensing Hydrologic Stations to achieve long sequence flow from high-precision remote sensing images. With remote sensing images downscaled from 30 m to 10 m resolution and multisource data integration of optical satellites, radar satellites and UAVs, river flow can be calculated with high accuracy and composited according to the rationality of flow data. In this case, different information can be obtained from multi-source remote sensing images, which indeed raises the utilization efficiency of these data sources. The method was tested on the Yunnan-Guizhou Plateau in southwest China and the accuracy of the results was verified. The results showed that (1) after this method was applied, the spatial resolution of Landsat-8 was significantly improved by linear downscaling and the water continuity index was increased by 7.1%. The average river length percentage in the channel reached 90.7%, an increase of 40%; (2) with the help of multisource data, the outliers of the flow series were greatly reduced, and the timespan was increased from approximately 10 data points in 5 years to more than 54 monthly flow data points, ensuring the temporal continuity of the calculated flow; (3) after downscaling, the images with better water continuity showed more stable river width, which reduced the RMSE of the Landsat-8 estimated flow from 13.5 m 3 /s to 5.4 m 3 /s, only 2/5 of the original. The MPE was reduced from 216% to 44%, specifically from 10.9 m 3 /s to 4.22 m 3 /s. The RMSE and MPE of the synthetic monthly sequence were also reduced, from 11.7 m 3 /s to 7.4 m 3 /s and from 10.0 m 3 /s to 5.99 m 3 /s respectively, after this method was applied. Our method is a complete set of mechanisms for image acquisition, downscaling, flow estimation and composited adjustment. It can be applied universally for flow estimation of small-and medium-scale rivers. Our method can provide a means of flow estimation for areas where the scope of remote sensing methods has previously been limited. In addition, this method can directly obtain continuous and high-precision water area inundation maps without the interference of weather and time, which will be of great value for water resource management and utilization.