Secondary Precipitation Estimate Merging Using Machine Learning: Development and Evaluation over Krishna River Basin, India

: The study proposes Secondary Precipitation Estimate Merging using Machine Learning (SPEM2L) algorithms for merging multiple global precipitation datasets to improve the spatiotemporal rainfall characterization. SPEM2L is applied over the Krishna River Basin (KRB), India for 34 years spanning from 1985 to 2018, using daily measurements from three Secondary Precipitation Products (SPPs). Sixteen Machine Learning Algorithms (MLAs) were applied on three SPPs under four combinations to integrate and test the performance of MLAs for accurately representing the rainfall patterns. The individual SPPs and the integrated products were validated against a gauge-based gridded dataset provided by the Indian Meteorological Department. The validation was applied at di ﬀ erent temporal scales and various climatic zones by employing continuous and categorical statistics. Multilayer Perceptron Neural Network with Bayesian Regularization (NBR) algorithm employing three SPPs integration outperformed all other Machine Learning Models (MLMs) and two dataset integration combinations. The merged NBR product exhibited improvements in terms of continuous and categorical statistics at all temporal scales as well as in all climatic zones. Our results indicate that the SPEM2L procedure could be successfully used in any other region or basin that has a poor gauging network or where a single precipitation product performance is ine ﬀ ective.


Introduction
Measurement of precipitation, an intrinsic component of the water cycle is highly essential for hydrological, climatological, agricultural, and environmental modeling purposes such as rainfall-runoff modeling, extreme event detection, weather prediction, hydraulic structure design, and water quality and quantity planning and management [1]. Precipitation can be derived from gauge stations, ground-based weather radars, spaceborne weather radars, and satellite radiometers (infrared and microwave) at various spatial and temporal scales [2]. The gauge-measured precipitation estimates that are considered accurate may lack representativity of surface rainfall over certain regions due to inadequate coverage over hilly terrains and unavailability over oceans. The rain gauge measurements even suffer from specific errors related to wind speed, evaporation and wetting processes [3].
Weather radars, on the other hand, employ active remote sensing sensors for measuring precipitation over large areas. The complex processing involved for converting the data obtained from these active remote sensing sensors induces uncertainty related to signal attenuation, range degradation, residual clutter, and variation of measured reflectivity over various spatial ranges [4]. Though the satellite radiometer estimates are not as accurate as gauge based or weather radar precipitation estimates, they provide global scale high-resolution precipitation estimates at various spatial and temporal scales. These satellite radiometer precipitation estimates are highly useful in developing countries that lack dense rain gauge network and radar infrastructure [2]. These satellite precipitation estimates (purely derived from satellite sensors, e.g., SoilMoisture2RAIN-Climate Change Initiative (SM2RAIN-CCI dataset)) and secondary satellite precipitation estimates (either satellite data adjusted with gauge data or reanalysis datasets) have been widely applied over a range of studies related to hydrological modeling [5][6][7][8][9], rainfall erosivity mapping [10][11][12], rainfall characterization [13,14], ecology [15][16][17], forest fire mapping [18,19], natural disasters [20][21][22], and soil moisture [23][24][25].
Recently, apart from these statistical merging techniques, very few studies employed machine learning models (MLMs) for fusing different satellite datasets related to precipitation, soil moisture, evapotranspiration (ET) and topography to improve the spatial variability of rainfall [36][37][38][39]. Some of these studies implemented one MLA to merge different datasets, or the precipitation datasets were merged with other hydrological variables (ET and soil moisture) and topographic datasets (slope and elevation). To date, no study tested more than one MLA for merging different rainfall products in any basin. For improving the accuracy of secondary precipitation products (SPPs) while preserving their spatial characterization of rainfall fields, it is possible to merge multiple SPPs, and doing so enables integrating the advantages of each sensor/algorithm, while partially overcoming their drawbacks. As per our knowledge, no study was found related to the merging of secondary rainfall products alone (not including soil moisture, ET, and topographic variables) implementing MLAs. Merging SPPs alone might overcome the underestimations (usually seen in most of the SPPs that are derived from IR sensors due to fixed brightness thresholds) and biases (errors) in individual datasets resulting in an improved merged dataset that consists of limited bias and underestimation. Further, as precipitation forms base for many hydrological, meteorological, and environmental applications, developing an accurate merged precipitation product by integrating data obtained from different earth observation satellite sensors employing MLAs is essential. The developed merged product can further be used for predicting disasters such as droughts and floods accurately to reduce risk and to develop resilience. As United Nation's Sustainable Development Goal (SDG) 15 concentrates on monitoring and identifying challenges and concerns in monitoring and efficiently maintaining environmental resources [40], the methodology can be implemented in any region, basin, or country for developing real-time integration systems and decision support systems for planning and developing sustainable management strategies.
Based on the research gaps discussed above, the current study was formulated by employing sixteen MLAs nested under six machine learning models (MLMs) to merge three different SPPs under four combinations to test the improvement in the spatial variability of rainfall over the Krishna River Basin (KRB), India. The current objectives of the present study are (1) testing the individual performance of Climate Hazards Group Infrared Precipitation with Station data (CHIRPS), Precipitation Estimation from Remotely Sensed Information using Artificial Neural Networks-Climate Data Record (PERSIANN-CDR), and recently released ECMWF (European Centre for Medium-Range Weather The major tributaries of KRB are Tungabhadra, Malaprabha, Koyna, and Bhima rivers [42]. Dominant soils found in this region are black soils, mixed soils, saline and alkaline soils, alluvium, and red soils. Much of the basin falls under a semi-arid climatic zone with massive occupancy of agricultural lands (76%). The majority of rainfall (90%) occurs during monsoon (June-October), with an average annual rainfall of 780 mm [43]. The maximum temperature varies between 20° and 42 °C, and the minimum temperature varies between 8° and 30 °C. The entire basin can be classified into The major tributaries of KRB are Tungabhadra, Malaprabha, Koyna, and Bhima rivers [42]. Dominant soils found in this region are black soils, mixed soils, saline and alkaline soils, alluvium, and red soils. Much of the basin falls under a semi-arid climatic zone with massive occupancy of agricultural lands (76%). The majority of rainfall (90%) occurs during monsoon (June-October), with an average annual rainfall of 780 mm [43]. The maximum temperature varies between 20 • and 42 • C, and the minimum temperature varies between 8 • and 30 • C. The entire basin can be classified into three climatic zones as the Tropical Monsoon rainforest (Am), Tropical savannah climate with Winter Dry (Aw), and Semi-Arid steppe climate (BS) classification as per Koppen climate classification system. More details on the classification criteria implemented to segregate into different climatic zones in the Koppen system can be found in Chen and Chen (2013) and Kottek et al. (2006) [44,45].

Datasets Used
Four gridded precipitation products having 0.25 • × 0.25 • spatial resolution with a temporal resolution spanning from 1985 to 2018 were implemented in the current study. The four gridded SPPs employed in the present study are described below. The main reason for selecting these datasets is due to their long term and real-time availability, finer spatial resolution and their consistent performance that was reported in various studies [46][47][48][49][50].
The recently released ERA-5 dataset was developed from the ECMWF atmospheric reanalysis product. The dataset is available at a daily time step for a global scale with a resolution of 0.25 • × 0.25 • . The dataset was developed employing reanalysis procedures by implementing data from various model simulations across the world using different laws of physics. More details regarding the dataset development can be found in Copernicus Climate Change Service, 2017 [51]. The data was downloaded from https://cds.climate.copernicus.eu/.

Indian Meteorological Department (IMD)
IMD has developed a high-resolution daily rainfall gridded dataset of 0.25 • × 0.25 • for the entire Indian subcontinent. The Shepard Interpolation technique was implemented for interpolating all the precipitation values obtained from individual gauge stations. IMD used~1803 gauge rainfall stations data for estimating 24-h accumulated rainfall. More details regarding IMD gridded data can be found in Pai et al. (2014) [52], and the data can be accessed from http://www.imdpune.gov.in/.

CHIRPS
Satellite and gauge-based rainfall estimates are provided by Climate Hazards Group Infrared Precipitation (CHIRP) product that has finer spatial resolution and longer periods. The CHIRPS (CHIRP blended with station data) datasets are developed by the University of California and USGS Geologic Survey that has a finer spatial resolution of 0.05 • covering 50 • N-50 • S spanning from 1981. Different products like Climate Hazards Group Precipitation Climatology (CHPClim) having 0.05 • resolution at monthly time scale estimated from station data, thermal infrared-based satellite rainfall estimates, Tropical Rainfall Measuring Mission (TRMM 3B42 v7) rainfall data, National Oceanic and Atmospheric Administration (NOAA) climate forecast system atmospheric rainfall data, and station-based precipitation data are blended and processed to produce CHIRPS rainfall data with high spatial and temporal resolution. The product so formed is of the penta-daily time step, which is then aggregated to a monthly time step. The present study uses CHIRPS data with 0.25 • spatial resolution spanning for the years 1985 to 2018. More details regarding IMD gridded data can be found in Funk et al. (2015) [53], and the data can be accessed from the CHG portal or the mentioned FTP site (ftp://ftp.chg.ucsb.edu/pub/org/chg/products/).

Precipitation Estimation from Remotely Sensed Information Using Artificial Neural Networks-Climate Data Record (PERSIANN-CDR)
The need for long-term rainfall records is addressed by the development of the PERSIANN-CDR product. PERSIANN-CDR has a spatial resolution of 0.25 • covering 60 • S to 60 • N having a span from January 1983 to December 2018. Global precipitation estimates provided by the PERSIANN algorithm uses both Infrared satellite and passive microwave data from GEO and LEO satellites. The 3-h Remote Sens. 2020, 12, 3013 5 of 23 rainfall estimates are generated by applying the PERSIANN algorithm to GridSat B1 Infrared satellite observations. The PERSIANN algorithm finds out the relationship between microwave and infrared data based on a neural network approach, and this relation is applied to Infrared (IR) data to generate 3-h rainfall estimates. The 3-h precipitation estimates are then combined with the GPCP product to eradicate the biases in the product that produces an adjusted PERSIANN-B1 rainfall estimate. The so-formed 3-h precipitation estimates after bias correction are then accumulated to form daily PERSIANN CDR Product of spatial resolution 0.25 • . More details regarding PERSIANN-CDR data can be found in Ashouri et al. (2015) [54], and the data can be accessed from http://chrsdata.eng.uci.edu/ portal.

Implementation of Machine Learning Techniques
In the current study, Secondary Precipitation Estimate Merging using Machine Learning (SPEM2L) was carried out to find out the optimum MLA that can be employed for merging precipitation datasets to generate improved spatial precipitation values over a basin. Sixteen MLAs broadly classified under six MLMs represented in Figure 2 were implemented in the current study. All the MLMs were trained using MATLAB 2018a, on an Intel Core i7, 3.07 GHz machine with 16 GB RAM. As this is a preliminary study to find the best MLA out of many existing algorithms to merge different precipitation datasets effectively, all the MLAs are trained with default hyperparameters. The hyperparameters are used to control the learning process of MLAs that are required to be specified before starting the learning process in a training phase. The hyperparameters are used to build relations or to fit an equation between input and output variables. For example, in the case of neural networks, the hyperparameters are the number of hidden units, activation function, number of neurons in each layer, etc. The hyperparameters employed in the current study are represented in Table 1. Various MLMs represented in Figure 2 were developed and evaluated in the selected study area.
Remote Sens. 2020, 12, x FOR PEER REVIEW 5 of 23 infrared data based on a neural network approach, and this relation is applied to Infrared (IR) data to generate 3-hr rainfall estimates. The 3-hr precipitation estimates are then combined with the GPCP product to eradicate the biases in the product that produces an adjusted PERSIANN-B1 rainfall estimate. The so-formed 3-hr precipitation estimates after bias correction are then accumulated to form daily PERSIANN CDR Product of spatial resolution 0.25°. More details regarding PERSIANN-CDR data can be found in Ashouri et al. (2015) [54], and the data can be accessed from http://chrsdata.eng.uci.edu/ portal.

Implementation of Machine Learning Techniques
In the current study, Secondary Precipitation Estimate Merging using Machine Learning (SPEM2L) was carried out to find out the optimum MLA that can be employed for merging precipitation datasets to generate improved spatial precipitation values over a basin. Sixteen MLAs broadly classified under six MLMs represented in Figure 2 were implemented in the current study. All the MLMs were trained using MATLAB 2018a, on an Intel Core i7, 3.07GHz machine with 16 GB RAM. As this is a preliminary study to find the best MLA out of many existing algorithms to merge different precipitation datasets effectively, all the MLAs are trained with default hyperparameters. The hyperparameters are used to control the learning process of MLAs that are required to be specified before starting the learning process in a training phase. The hyperparameters are used to build relations or to fit an equation between input and output variables. For example, in the case of neural networks, the hyperparameters are the number of hidden units, activation function, number of neurons in each layer, etc. The hyperparameters employed in the current study are represented in Table 1. Various MLMs represented in Figure 2 were developed and evaluated in the selected study area.

Linear Regression based Models
Simplified Linear Regression (LR), Linear Regression algorithm with Interactions term (LR-I), Stepwise Linear Regression algorithm (LR-S), and Linear Regression algorithm with Robust fitting (LR-R) were categorized and used under LR based models. A "fitlm" function with different options in MATLAB was used for LR, LR-I, and LR-R algorithms. A simple linear regression algorithm with stepwise regression was implemented for the LR-S algorithm. The default ordinary least squares method was employed in the case of LR to fit the algorithm between secondary precipitation datasets Remote Sens. 2020, 12, 3013 6 of 23 (CHIRPS + ERA-5) and gauge-based gridded data (IMD). A "bisquare" robust fitting weighting function was implemented with default tuning constant for the LR-R algorithm to fit between satellite and gauge based gridded datasets. Along with intercept and linear terms, an interaction term is further added to fit the LR-I MLA. Table 1. Hyperparameter values for all the 16 machine learning algorithms used in the study.

1.
Linear Regression with Robust Fitting Bisquare weight function with a tuning constant of 4.685 2.

k-Nearest Neighbour
Nearest neighbour search method-kdtree Distance calculation method-Euclidean Distance weighting function-Equal Method of breaking ties-smallest Bucket Size-50 (maximum number of data points in the leaf node of the kd-tree) 3.

Ensemble Trees (For both Bagged and Boosted methods)
Number of learning cycles-100 Learn Rate-1 4.

Neural Network (Levenberg Marquardt Optimization & Bayesian Regularization) Maximum number of epochs-1000
Minimum gradient-1 × 10 −7 Performance Function-RMSE Number of hidden layers-1 Activation functions-Hyperbolic tangent (tansig) for hidden layer and linear (purelin) transfer function for output layer Number of neurons in a hidden layer-2 (two datasets integration) and 3 (three datasets integration) Type of connection-dense connection Architecture-Feedforward neural network 5.

Support Vector Regression (Fine, Medium and Coarse) Box constraint (Maximum limit for Alpha coefficients) is Interquartile range of response variable/1.349 Gaussian or Radial Basis Function (RBF) kernel Optimization routine is Sequential Minimal Optimization
Maximum number of optimization iterations-1 x 10 6 . Kernel Scale for Fine, Medium, and Coarse: 0.61, 2.4, and 9.8, respectively.

Support Vector Machine Regression Models
SVMs are supervised learning techniques developed for classification problems; however, they can be implemented for regression as well. SVM has a significant advantage to avoid trapping in a local minimum and network overfitting. In SVM, a hyperplane is constructed in a higher-dimensional space to classify the data. A nonlinear function in SVM is implemented to map the input data to the higher dimensional space. Based on the classification technique of SVM, Support Vector Regression (SVR) is proposed. The error between the estimate and actual value is calculated using cost function in SVR, whereas a threshold value (є) is used to specify the distance between the boundary lines. Four SVR algorithms-SVR with linear kernel regression algorithm (SL), SVR with a Gaussian kernel regression function having a kernel scale of 0.61 for SVR with Fine Gaussian kernel (SF), 2.4 for SVR with medium Gaussian kernel (SM), and 9.8 for SVR with coarse Gaussian kernel (SC)-were implemented in the study. The function "fitrsvm" in MATLAB was used for SL, SF, SM, and SC algorithms [55].

Regression Tree Models
In the case of regression trees, each dataset has undergone recursive binary portioning at the nodes, and the relationship between the dependant response and predictor variables was explored. These regression tree-based algorithms, namely, Fine Tree (FT), Medium Tree (MT), and Coarse Tree (CT), were implemented using "fitrtree" function in MATLAB with minimum leaf sizes of 4, 12, and 36, respectively [56].

Neural Network Models
Neural Network-based Levenberg-Marquardt algorithm (NLM) and Bayesian Regularization algorithm (NBR) were constructed for minimizing the errors by employing the Multilayer Perceptron Neural network (MLPNN) architecture with Levenberg Marquardt [59,60] and Bayesian Regularization functions [61]. From the entire input data, 60% was used for training the networks, and the remaining 40% was used for testing and validation purposes. The MLPNN architecture consists of single input, hidden and output layers. The hidden layer consists of two/three neurons corresponding to two/three inputs of secondary precipitation datasets [62].

K-Nearest Neighbour Models
The "fitcknn" function in MATLAB with dependent options is used in the current study. The regression fit between SPPs and IMD gridded data was carried out by employing a single neighbor and Euclidean distance in the current study [63,64].

SPEM2L Procedure
The selected SPPs were divided into two groups: a training set (to train the MLMs) and a validation set (to evaluate the performance of the integrated product). In the current study, 34 years of precipitation data spanning from 1985 to 2018 was considered for merging precipitation datasets using MLAs. Out of these 34 years, 22 years (1985-2006) of data were considered for training the MLAs, whereas the other 12 years (2007-2018) of rainfall data were employed for testing and validation purposes. The continuous rainfall data was considered for training and testing to ensure that the minimum and maximum precipitation values that occur over summer and monsoon seasons are equally covered during training, testing, and validation phases. All the precipitation datasets are resampled to the same spatial (0.25 • ) and temporal resolution (daily) and to the same grid location as of IMD to ensure similar raster geometry. The SPPs were resampled by employing a cubic spline interpolation algorithm in MATLAB. During the training phase, all the MLMs were trained on the training subset implementing supervised regression between the SPPs and the IMD gridded data. As three SPPs are selected as an input, four combinations are developed to merge different SPPs and to evaluate against IMD gridded data: • Combination 1: CHIRPS and ERA-5 precipitation datasets were merged employing sixteen MLAs and evaluated against IMD gridded data. To carry out these combinations, first, the precipitation values at each grid cell are extracted for the entire period . Second, for each day, an MLA was trained for the years 1985-2006 using the selected SPPs in each combination as the dependent variable, and the precipitation values from IMD data as predictors. Third, the trained MLA was used to predict the daily precipitation values for each pixel of the study area for the remaining years (2007)(2008)(2009)(2010)(2011)(2012)(2013)(2014)(2015)(2016)(2017)(2018). This process was repeated for each day for the years 1985-2018, employing all sixteen algorithms.

Performance Evaluation
Statistical analysis was performed to compare each SPP and the integrated products against IMD gridded data by employing categorical and continuous statistics that communicates the detection capabilities and error characteristics of SPPs. Probability of Detection (POD), False Alarm Ratio (FAR) categorized under categorical statistics, Correlation Coefficient (R), and Root Mean Square Error (RMSE) grouped under continuous statistics were employed in the current study for comparison. The categorical metrics convey the precipitation detection potentiality of SPPs. In contrast, the continuous metrics disclose the performance of SPP in estimating the amount of precipitation accurately when compared with the standard dataset (IMD). The ratio of hits (when SPP accurately detected the rainfall event concerning standard dataset) to the actual number of precipitation events occurred according to the base dataset is represented by POD. The ratio of misses (when SPP detects the rainfall event during the absence of rainfall in the standard dataset) to the actual number of rainfall events that are not diagnosed by the reference dataset is represented by FAR. A threshold of 1 mm/day was implemented for computing POD and FAR in the study. The categorical metrics were computed for the entire time series and it was segregated for different rainfall regimes. The precipitation intensities were partitioned into multiple rainfall regimes as low, medium, and high based on the criteria mentioned in Table 2. Table 2. Partition of overall precipitation time series into low, medium, and high precipitation time series (PCP is precipitation, µ is mean of rainfall, and σ is the standard deviation of precipitation).

Rainfall Regime Criterion
Low PCP < µ Medium PCP ≥ µ and PCP ≤ µ + 2σ High PCP > µ + 2σ The degree of relevance in precipitation between SPP and the gauge-based dataset is represented by R, RMSE indexes, or the averaged error magnitude between the SPP and the standard dataset. Higher R and lower RMSE represents greater accuracy of SPP in terms of detecting precipitation concerning base or standard datasets. More details regarding the metric and formulas can be found in

Spatial Pattern Assessment of the Rainfall Products
The spatial maps of average annual precipitation computed employing 34 years of precipitation data from IMD, CHIRPS, PERSIANN-CDR, and ERA-5 datasets for the Krishna River Basin in India are represented in Figure 3. From Figure 3, it can be observed that the average annual rainfall recorded by standard dataset (IMD) varies from 300 to 3000 mm throughout the basin. From Figures 1 and 3, it can be noticed that high rainfall patterns were concentrated over the hilly regions of the basin that are having an elevation higher than 1000 m above mean sea level. Low rainfall pattern was recorded in the central part of the basin that falls under the semi-arid climatic zone (Figures 1 and 3). A medium amount of rainfall ranging from 1000 to 1500 mm was concentrated at the mouth portion of the basin. Figure 3 illustrates that the western part of the basin classified under the tropical climatic zone recorded the highest amount of rainfall in all datasets due to the orographic effect of precipitation that occurs Remote Sens. 2020, 12, 3013 9 of 23 due to the occupancy of the Western Ghats mountains in the upper part of the basin. PERSIANN-CDR underestimated the amount of rainfall throughout the basin when compared to IMD with the average annual rainfall ranging from 800 to 1400 mm. The ERA-5 dataset contributed to an equal amount of precipitation (approximately 30-3000 mm) concerning IMD over the entire watershed. The CHIRPS dataset reported an underestimation with lower precipitation amounts in the central portion of the basin compared to the IMD dataset ( Figure 3).          (Figures 4 and 5). The comparable performance was achieved by LR, LR-I, LR-S, and NLM next to NBR with a median R-value of 0.33 and RMSE of 11.77 mm/day. From the individual comparison results of three SPPs with IMD, ERA-5 proved better when compared to CHIRPS and PERSIANN-CDR with lower error (RMSE) and higher correlation (R). The combinations that involved ERA-5 (CHIRPS + ERA-5 or PERSIANN + ERA-5) exhibited better performance than the combination that is having the other two datasets (CHIRPS + PERSIANN-CDR).

Evaluation of Machine Learning Models Performance at Daily Time Step
From the results of four combinations, it can be clearly stated that NBR outperformed from all MLMs in all combinations when either two or three datasets are merged. NLM and LR resulted in a comparable performance next to NBR in all combinations (Figures 4 and 5). The KNN algorithm proved ineffective for combining different rainfall datasets in all sets of combinations. These results suggest that neural networks have high capability in merging different rainfall datasets. The simplest LR even exhibited comparable performance with NBR, suggesting that LR can be applied if complex trees or neural networks cannot be modelled and applied. From these results, it can be concluded that the SPEM2L procedure was able to improve the spatiotemporal representation of precipitation by integrating multiple rainfall datasets. Our results exhibited that the blending of multiple precipitation estimates can improve the spatiotemporal characterization of rainfall, which is consistent with the results obtained by Manz et al. (2016) and Verdin et al. (2016) [68,69]. From the overall results, the NBR algorithm from combination 4 exhibited better performance when compared to the other 15 MLAs and three combinations. As NBR from combination 4 (hereafter called as NBRC-4) proved effective, further analysis in the manuscript was performed and discussed with NBRC-4 results.

Categorical Metric Assessment
As the NBRC-4 algorithm performed better compared to individual datasets (ERA-5, PERSIANN-CDR, and CHIRPS) and other model combinations in continuous statistic assessment, the categorical statistical metric evaluation is discussed with the results obtained from NBRC-4 algorithm. Categorical statistics, i.e., POD and FAR, are computed based on the formulas mentioned in Table 2 and are segregated into low (0-5 mm), medium (5-25 mm), and high (>25 mm) rainfall Remote Sens. 2020, 12, 3013 14 of 23 classes. The pattern of POD and FAR for the NBRC-4 algorithm dataset computed against IMD is shown in Figure 6. The NBRC-4 algorithm has a high POD value (median of 0.572) in low rainfall criteria indicating better rainfall event detection at a daily time step. For medium and high rainfall intensities, the value of POD followed a declining trend (median of 0.09 and 0.009), indicating that the rainfall event detection capability decreases with the increase of rainfall intensity, as shown in Figure 6.
Remote Sens. 2020, 12, x FOR PEER REVIEW 14 of 23 in Table 2 and are segregated into low (0-5 mm), medium (5-25 mm), and high (>25 mm) rainfall classes. The pattern of POD and FAR for the NBRC-4 algorithm dataset computed against IMD is shown in Figure 6. The NBRC-4 algorithm has a high POD value (median of 0.572) in low rainfall criteria indicating better rainfall event detection at a daily time step. For medium and high rainfall intensities, the value of POD followed a declining trend (median of 0.09 and 0.009), indicating that the rainfall event detection capability decreases with the increase of rainfall intensity, as shown in Figure 6. FAR, which should have a value closer to 0, indicated moderate rainfall detection capability at low rainfall events with a median value of 0.38. The NBRC-4 algorithm detected more false rainfall events in case of medium rainfall resulting in a median FAR value of 0.48. In the high rainfall intensity category, NBRC-4 performed better, detecting fewer false rainfall events (median value of 0.29) when compared with the IMD dataset ( Figure 6). The POD values with 1mm/day rainfall threshold exhibited that the overall rainfall detection skill of all datasets decreases as precipitation intensity increases, and vice versa for FAR metrics. These categorical assessment results reveal that all the implemented datasets cannot detect the magnitude of extreme precipitation events accurately. Similar results were reported by Ahmed et al. (2019) and Pour et al. (2014) [70,71], where they concluded that high magnitude or extreme rainfall occurs at a micro-level and therefore, the SPPs may not have the potentiality to detect those rainfall events accurately at the point level. In supporting this statement, Schneider et al. (2017) [72] commented that the SPPs could detect these extreme rainfall events when the gridded dataset is developed with forcing a more significant or dense number of observed station data over an area which we currently lack all around the world. From the overall categorical skill score test results, it can be observed that the NBRC-4 algorithm performed better in detecting low rainfall events compared to medium and high-intensity precipitation events.

Temporal Assessment of NBRC-4 Algorithm Performance
The performance of the NBRC-4 algorithm is assessed at a 3-day and monthly time step by applying continuous statistic metrics (R and RMSE). The performance of individual datasets (ERA-5, PERSIANN-CDR, and CHIRPS) at a 3-day time step is similar to that of a daily time step that has a median R-value of 0.21, 0.28, and 0.19 and median RMSE value of 15.38, 13.2, and 15.64 mm/3-day respectively (Figure 7). The NBRC-4 algorithm proved best even at a 3-day time step with a median R of 0.32 and RMSE of 11.28 mm. At a monthly time frame, individual datasets (ERA-5, PERSIANN-CDR, FAR, which should have a value closer to 0, indicated moderate rainfall detection capability at low rainfall events with a median value of 0.38. The NBRC-4 algorithm detected more false rainfall events in case of medium rainfall resulting in a median FAR value of 0.48. In the high rainfall intensity category, NBRC-4 performed better, detecting fewer false rainfall events (median value of 0.29) when compared with the IMD dataset ( Figure 6). The POD values with 1 mm/day rainfall threshold exhibited that the overall rainfall detection skill of all datasets decreases as precipitation intensity increases, and vice versa for FAR metrics. These categorical assessment results reveal that all the implemented datasets cannot detect the magnitude of extreme precipitation events accurately. Similar results were reported by Ahmed et al. (2019) and Pour et al. (2014) [70,71], where they concluded that high magnitude or extreme rainfall occurs at a micro-level and therefore, the SPPs may not have the potentiality to detect those rainfall events accurately at the point level. In supporting this statement, Schneider et al. (2017) [72] commented that the SPPs could detect these extreme rainfall events when the gridded dataset is developed with forcing a more significant or dense number of observed station data over an area which we currently lack all around the world. From the overall categorical skill score test results, it can be observed that the NBRC-4 algorithm performed better in detecting low rainfall events compared to medium and high-intensity precipitation events.

Temporal Assessment of NBRC-4 Algorithm Performance
The performance of the NBRC-4 algorithm is assessed at a 3-day and monthly time step by applying continuous statistic metrics (R and RMSE). The performance of individual datasets (ERA-5, PERSIANN-CDR, and CHIRPS) at a 3-day time step is similar to that of a daily time step that has a median R-value of 0.21, 0.28, and 0.19 and median RMSE value of 15.38, 13.2, and 15.64 mm/3-day respectively (Figure 7). The NBRC-4 algorithm proved best even at a 3-day time step with a median R of 0.32 and RMSE of 11.28 mm. At a monthly time frame, individual datasets (ERA-5, PERSIANN-CDR, and CHIRPS) exhibited a high correlation with a median R-value of 0.76, 0.75, and 0.75 and higher RMSE values of 53.36, 55.76, and 70.7 mm/month, respectively. However, the NBRC-4 algorithm outperformed even at a monthly time step with a median R and RMSE values of 0.84 and 40.5 mm/month, respectively, as shown in Figure 7. The evaluated merged precipitation products exhibited higher performances at the monthly time scales in comparison to shorter temporal scales (daily), similar to the results reported by Jiang et al. (2012) and Zambrano-Bigiarini et al. (2017) [73,74]. This represents that despite systematic, random, and detection errors present in precipitation datasets at the daily time step, they are still able to describe precipitation patterns when aggregated at longer temporal scales (monthly).  [73,74]. This represents that despite systematic, random, and detection errors present in precipitation datasets at the daily time step, they are still able to describe precipitation patterns when aggregated at longer temporal scales (monthly).

Performance Assessment of NBRC-4 Algorithm in Various Climatic Zones
The Krishna River Basin (KRB) employed in the current study comprises three climatic zones, i.e., Tropical Monsoon rainforest (Am), Tropical savannah climate with Winter Dry (Aw), and Semi-Arid Steppe climate (BS), as shown in Figure 1.
The performance of the NBRC-4 algorithm was assessed in the three climatic zones using continuous statistics. From

Performance Assessment of NBRC-4 Algorithm in Various Climatic Zones
The Krishna River Basin (KRB) employed in the current study comprises three climatic zones, i.e., Tropical Monsoon rainforest (Am), Tropical savannah climate with Winter Dry (Aw), and Semi-Arid Steppe climate (BS), as shown in Figure 1.
The performance of the NBRC-4 algorithm was assessed in the three climatic zones using continuous statistics. From   [77][78][79]while using CHIRPS dataset and the reason for the weaker performance may be attributed to the Infrared algorithm (which computes data from Infrared region retrieved signals) that was used for CHPClim dataset developments which implement fixed brightness temperature thresholds to differentiate raining and non-raining clouds. The defined thresholds are usually too cold and the orographic precipitation occurring over the Western Ghats was  [77][78][79] while using CHIRPS dataset and the reason for the weaker performance may be attributed to the Infrared algorithm (which computes data from Infrared region retrieved signals) that was used for CHPClim dataset developments which implement fixed brightness temperature thresholds to differentiate raining and non-raining clouds. The defined thresholds are usually too cold and the orographic precipitation occurring over the Western Ghats was warm which may not produce much ice aloft resulting in an underestimation of rainfall [77,[80][81][82][83][84]. The pattern of the RMSE statistical metric is like that of R, where ERA-5 revealed better results than CHIRPS and PERSIANN-CDR. From the overall results of continuous statistics at a daily time step, ERA-5 yielded better results than CHIRPS and PERSIANN-CDR in all climatic regions with high R and low RMSE values. From the overall results related to the assessment of the NBRC-4 algorithm in contrasting climatic zones, it is observed that the NBRC-4 algorithm outperformed when compared to individual dataset statistical performance.
Although the SPEM2L procedure was only applied over the KRB, we are confident that this method could be successfully applied over other basins or regions that are ungauged or have a poor gauging network, due to its effective performance in a region with notable heterogeneity in topography and climate; and because it was able to improve the spatiotemporal characterization of precipitation. Future studies can implement hyperparameter tuning and optimization for obtaining more effective results when implementing the NBR algorithm for merging different datasets. Researchers can further implement and test this NBR algorithm in different climatic and topographic regions for merging different SPPs to frame a standard framework for generating accurate precipitation products that can be used for monitoring and predicting natural hazards. Apart from the employed SPPs in the current study, different precipitation products should be implemented and merged based on their availability to test the efficiency of NBR algorithm for merging different SPPs. The merged precipitation product can further be implemented in a hydrological model to test its efficiency in simulating water balance components accurately compared to individual precipitation dataset simulations.

Conclusions
Secondary precipitation products furnish an exceptional anomalous opportunity for various hydro-meteorological and environmental applications. Despite the continuous improvements of rainfall datasets, still, the datasets are prone to biases and mismatches. To overcome some of those biases, here we presented a novel SPEM2L framework that is capable of deriving improved precipitation estimates by merging secondary rainfall datasets alone that are publicly available for the entire globe with near-real-time precipitation estimates. Neural network-based algorithms (NLM and NBR) and linear regression-based algorithms (LR, LR-I, and LR-S) developed with the proposed method exhibited improved rainfall characterization at all temporal scales and contrasting climatic zones compared to the individual precipitation datasets. Finally, the NBR algorithm with three dataset integration performed better than other combinations and other 15 MLA's. Some of the key findings from this applied methodology are mentioned below.

1.
ERA-5 yielded better statistical results than CHIRPS and PERSIANN-CDR in all climatic zones with high correlation and low magnitude of errors (RMSE).

2.
The POD and FAR statistics exhibited that the overall rainfall detection skill of all datasets decreases as precipitation intensity increases.

3.
The individual precipitation products and the MLA's exhibited superior performance at longer timescales (monthly) than at shorter spans (daily). 4.
NBRC-4, which employs three-dataset merging, performed better than all other combinations that involve two rainfall dataset integration. 5.
NBRC-4 algorithm outperformed than all the evaluated rainfall products at different temporal resolutions (daily, 3-day, and monthly) and climatic zones (Am, Aw, and BS). 6.
The amplification of RMSE can be observed from arid to tropical climatic zones in all the tested precipitation datasets. 7.
SPEM2L procedure can be applied at different temporal scales (i.e., daily, 3-day, and monthly) to acquire an improved spatiotemporal rainfall characterization.
SPEM2L was proposed and tested in a heterogeneous catchment (KRB) that has varying climate and topography with publicly available global datasets. The proposed methodology improved the spatiotemporal characterization of rainfall, and hence the authors are confident that this framework can be applied for any region or basin that has a poor gauging network or where a single precipitation product is testing ineffective. The current research gives an insight to the researchers to select suitable MLAs for integrating multiple precipitation datasets along with the best dataset that can be employed in any region for hydrological and climatological applications. Funding: This research received no external funding and the APC was funded by MDPI.

Acknowledgments:
The authors are very thankful to all the organizations and institutes that developed precipitation datasets and made them available for free.

Conflicts of Interest:
The authors declare no conflict of interest.