remote sensing Article Assessing Near Real-Time Satellite Precipitation Products for Flood Simulations at Sub-Daily Scales in a Sparsely Gauged Watershed in Peruvian Andes Harold Llauca 1,* , Waldo Lavado-Casimiro 1 , Karen León 1 , Juan Jimenez 1, Kevin Traverso 1 and Pedro Rau 2 1 Servicio Nacional de Meteorología e Hidrología del Perú (SENAMHI), Lima 15072, Peru; wlavado@senamhi.gob.pe (W.L.-C.); kleon@senamhi.gob.pe (K.L.); jjimenez@senamhi.gob.pe (J.J.); arnold.traverso@gmail.com (K.T.) 2 Departamento de Ingeniería Ambiental, Centro de Investigación y Tecnología del Agua (CITA), Universidad de Ingeniería y Tecnología (UTEC), Lima 15063, Peru; prau@utec.edu.pe * Correspondence: hllauca@senamhi.gob.pe; Tel.: +511-954301911 Abstract: This study investigates the applicability of Satellite Precipitation Products (SPPs) in near real-time for the simulation of sub-daily runoff in the Vilcanota River basin, located in the southeastern Andes of Peru. The data from rain gauge stations are used to evaluate the quality of Integrated Multi-satellite Retrievals for GPM–Early (IMERG-E), Global Satellite Mapping of Precipitation– Near Real-Time (GSMaP-NRT), Climate Prediction Center Morphing Method (CMORPH), and HydroEstimator (HE) at the pixel-station level; and these SPPs are used as meteorological inputs for  the hourly hydrological modeling. The GR4H model is calibrated with the hydrometric station of  the longest record, and model simulations are also verified at one station upstream and two stations Citation: Llauca, H.; downstream of the calibration point. Comparing the sub-daily precipitation data observed, the Lavado-Casimiro, W.; León, K.; results show that the IMERG-E product generally presents higher quality, followed by GSMaP-NRT, Jimenez, J.; Traverso, K.; Rau, P. CMORPH, and HE. Although the SPPs present positive and negative biases, ranging from mild Assessing Near Real-Time Satellite Precipitation Products for Flood to moderate, they do represent the diurnal and seasonal variability of the hourly precipitation in Simulations at Sub-Daily Scales in a the study area. In terms of the average of Kling-Gupta metric (KGE), the GR4H_GSMaP-NRT’ Sparsely Gauged Watershed in yielded the best representation of hourly discharges (0.686), followed by GR4H_IMERG-E’ (0.623), Peruvian Andes. Remote Sens. 2021, GR4H_Ensemble-Mean (0.617) and GR4H_CMORPH’ (0.606), and GR4H_HE’ (0.516). Finally, the 13, 826. https://doi.org/10.3390/ SPPs showed a high potential for monitoring floods in the Vilcanota basin in near real-time at the rs13040826 operational level. The results obtained in this research are very useful for implementing flood early warning systems in the Vilcanota basin and will allow the monitoring and short-term hydrological Academic Editor: Alberto Refice forecasting of floods by the Peruvian National Weather and Hydrological Service. Received: 27 January 2021 Keywords: GPM; GSMaP; GR4H; floods; Andes; Peru; Vilcanota Accepted: 19 February 2021 Published: 23 February 2021 Publisher’s Note: MDPI stays neutral 1. Introduction with regard to jurisdictional claims in published maps and institutional affil- Floods are one of the most frequent natural disasters in Peru [1], negatively impacting iations. the country at a socioeconomic level, affecting the physical and mental health of the popula- tion and causing severe damage to public and private property each year [2]. For example, heavy precipitation in January 2010 caused a great impact on the Andean region of Cusco. The Peruvian National Institute of Civil Defense (INDECI) reported that approximately 12,000 homes were affected. Total economic losses were estimated at USD 220 million, with Copyright: © 2021 by the authors. 55% corresponding to infrastructure, 35% to health, and 10% to agriculture and tourism [3]. Licensee MDPI, Basel, Switzerland. This article is an open access article Despite the efforts made in recent decades, floods remain a recurring problem and tend distributed under the terms and to be increasingly frequent in the context of a changing climate. Therefore, flood monitor- conditions of the Creative Commons ing and forecasting are key to reducing risk, improving decision-making and population Attribution (CC BY) license (https:// response capacity, and mitigating its negative impacts [4]. creativecommons.org/licenses/by/ Precipitation is the most common meteorological forcing of a hydrological model [5,6]. 4.0/). Precipitation data typically come from rain gauge stations; however, in mountainous and Remote Sens. 2021, 13, 826. https://doi.org/10.3390/rs13040826 https://www.mdpi.com/journal/remotesensing Remote Sens. 2021, 13, 826 2 of 18 difficult-to-reach areas such as the Peruvian Andean region [7], the density of stations is generally insufficient to adequately represent the high variability of precipitation. The costly implementation and maintenance of a rain gauge network with automated real-time data transmission makes it even more difficult to improve the rain gauge network in the Andes. In contrast, with the improvement of remote sensing techniques and precipita- tion reconstruction algorithms, satellite precipitation products (SPPs) have gained more attention in recent years for meteorological and hydrological applications due to their broad spatial coverage and fine spatiotemporal resolution [2,4,5,8–10], especially in large basins, with long periods of concentration and with little precipitation information in the Peruvian Andes [11–14]. However, the errors inherent to the SPPs must be evaluated with imperfect ground-station observations to improve short-range precipitation forecast and reduce uncertainty and error propagation in the forecast chain by combing meteorological and hydrological information [15]. Flood Early Warning Systems (FEWSs) are the main tools applied for flood risk as- sessment and damage reduction by flood forecasting using real-time information obtained from rain gauges and radar [16]. The use of SPPs in operational hydrological applications remains challenging mainly due to their spatiotemporal resolution and latency [10], espe- cially for monitoring floods caused by intense rains lasting less than 24 h [17]. In this sense, there are near real-time SPPs with short time intervals between data capture and gener- ation of the precipitation product. For example, the products Integrated Multi-satellite Retrievals for GPM-Early (IMERG-E), Global Satellite Mapping of Precipitation–Near Real Time (GSMaP-NRT), Climate Prediction Center Morphing Method (CMORPH), and Hy- droEstimator (HE) have latencies less than 24 h and fine spatial (5–10 km) and temporal (0.5–1 h) resolutions [18]. The literature contains studies demonstrating the hydrological application of SPPs for forecasting flows in data scarce basins using conceptual models such as HYMOD, GR4J, and MILc, and synthetic aperture radar images to validate the simulated flood extents [4,19]. Different SPP bias correction techniques and their impact on flood event estimation have also been evaluated [16,20,21], and wavelet-based machine learning models have been used for real-time flood forecasting from SPP data [22]. This study uses the conceptual GR4H hydrological model [23] to simulate sub-daily streamflows from SPP data for use in operational flood forecasting in near real-time. GR4H is a parsimonious and easy-to-implement model that has been widely used to evaluate sub-daily runoff in basins with different hydrological regimes worldwide [24–26]. Its semi-distributed application coupled with hydrological routing models also allows the sub-daily streamflows to be estimated for operational flood forecasting purposes. This study evaluates the suitability of applying near real-time SPPs for monitoring and forecasting floods in the Vilcanota River basin in the Cusco Region of Peru by way of (a) meteorological and (b) hydrological analyses. The first focuses on comparing SPP and rain gauge station data at the pixel-point level, while the second evaluates the hydrological performance in the bias-corrected SPP-forced simulated streamflows. For this, the current network of hydrometeorological stations in the study domain is used, and the conceptual GR4H model is used to simulate the hourly streamflows of the basin. The results obtained in this work will be very useful for implementing a FEWS in the Vilcanota River basin and will contribute to the generation of flood alerts and forecasts by the Peruvian National Weather and Hydrological Service. 2. Study Area and Data 2.1. Study Area The basin of the Vilcanota River delimited at Intihuatana km105 stream gauge station (Figure 1) lies in the Cusco Region, southeast of Peru, and it is part of the Andean region of the Amazon basin [27]. It has a drainage area of 9594 km2 and a mean annual discharge of 130 m3/s. Its relief is predominantly rugged, characteristic of the Cordillera Vilcanota, with a mean slope and elevation of 15.6◦ and 4176 masl, respectively. Annual precipitation ranges from 800 to 1000 mm. The Vilcanota basin is located in a transition zone with Remote Sens. 2021, 13, x FOR PEER REVIEW 3 of 18 Remote Sens. 2021, 13, 826 of 130 m3/s. Its relief is predominantly rugged, characteristic of the Cordillera Vil3coafn1o8ta, with a mean slope and elevation of 15.6° and 4176 masl, respectively. Annual precipitation ranges from 800 to 1000 mm. The Vilcanota basin is located in a transition zone with trop- troipcaicl,a sl,usbu-,b a-,nadn edxterxattrraotproicpailc callimclaimteas,t eas s,haosrht owret tw seetassoenas (oNno(vNemovbeemr btoe rMtoarMcha),r cahn)d, aan ldong a ldornyg pderryiopde druiordindg uthrien rgestht eofr tehset oyfeatrh e[2y8]e.a Trh[e2 8b]a. siTnh ceobnatisninuaclolyn tsiunfufearlsly frsoumff eflrosofdros mthat flosoedvsertehlayt hseavrmer epleyohpalerm’s hpeeaolpthle, ’psrhiveaatlet hp,rporpievrattye, panrodp reerstuy,lta innd larregsue letcionnloamrgiec elocsosneosm toi cthe loslsoecsalt ocotmhemluoncaitlyc oanmdm thuen iCtyusacnod rethgeioCnu [s2c9o]. rTehgeiorenfo[2re9,] .thTehree riesf aonre u, rthgeenret niseeadn tuor gimenptle- nemedentot aim mpolnemitoernint ga amnodn filtooordin fgoraencdasfltionogd syfosrteecma sitni nthges byastseinm inin qtuhaesib-aresainl tiinmqeu. Hasoi-wreeavler, timthee. Hlooww deevnesri,ttyh aenldow shdoernt sreitcyoardnidngsh toimrterse ocof rtdhein cgutrirmenets noefttwhoercku orrf ernatinn getawugoerk stoaftiroanins in gatuhgee bsatasitnio nmsaiknetsh tehbisa stainskm daikffeiscuthlti s[2ta8s].k Tdhieffirecfuolrte[,2 8th].e Tuhser eoffo rees,titmheatuesde hoofuersltyim satteldlite hopurrelycipsaittaetlliiotne preocdipuicttast iwonithp rfoindeu scptsawtioitthemfipneorsapl arteiosotelumtipoonr aisl rae sporloumtiiosninigs alpterronmatisivineg for altsehrnoratt-itveermfo rflsohoodr ts-itmerumlafltioonds. simulations. FigFuirgeur1e. 1L. oLcoactaiotinono foft htheeV Viliclcaannoottaa rriivveerr bbaassiinn aatt IINNTT ggaauuggee-s-statatitoionn. .PlPulvuivomiometreitcr iacnadn hdyhdyrodmroemtreict rnicetnweotwrko artk thate tshtuedy domain. study domain. In Itnh isthwiso rwko, arkr,e par erseepnrteasteinvetasttiuvde ystduodmya dinowmaasind ewfianse ddfeofirntehde efxotrr atchteio enx,tprarocctieosnsi,n pgr,o- ancdebssiainsgc,o arnredc tbioians ocfotrhreechtioounr loyf SthPeP sh,obuertlwy eSePnPtsh, ebreatwngeeenof tlhaet irtaundgees soof ulathtit1u2d.9eºs– s1o4u.8tºh, a1n2d.9º– lon1g4i.t8uº,d aensdw leosntg7i0tu.6dº–e7s 2w.7eºs. tA 7d0.d6iºt–io 7n2a.7llºy. ,A1d1dsuitbio-bnaaslliyn,s 1a1n sdurbi-vbearssinecst aionnds rwiveerre sdeectliinoenast wedere fordseelimnei-adtiesdtr ifbour tseedmhi-yddirsotrloibguicteadl m hoydderolilnoggiucsailn mg othdeelGinRg4 Husimngod thele( GFiRg4uHre m1)o, cdoenl s(iFdiegruinreg 1), thecolnosciadtieorningo fththee loecxaistitoinng ohf ythder oemxiesttrinicg shtaytdioronms ientritch estbataisoinns. iTnh tehed ebtaasiilns. oTfhteh edemtaaiilns of phtyhsei omgraainp hpihcycshiaorgarcatpehriisct icchsaorfacetaecrhisstiucbs- obfa seiancha rseusbh-obwasnini narTea sbhleow1.n in Table 1. 2.22..2G. rGouronudn-Bda-BseadseDd adtata HoHuorluyrldya dtaatwa ewreerceo lcloelclteecdtefdro fmrom15 1r5a irnaigna guaguegseta sttiaotniosnasn adnfdo uforuhry hdyrodmroemtreictrsicta sttiaotniosns frofmromth ethaeu atoumtoamteadtende ntwetowrkorokf othf ethPee PruervuiavniaNn aNtiaotnioanl aWl eWatehaetrhearn adnHd yHdyrodlroogloicgailcSael rSveircveice (SE(SNEANMAHMIH) If)o frotrh tehep epreiroidodfr formom1 1J aJannuuaarryy2 2001166 ttoo 31 March 2020.. TTaabbllee 22 ssuummmmaarirzizese sthe theseslelcetcetde dstastaiotinosn. sT.hTe hloeclaolc tailmtiem ise fiisvefi vheouhrosu lresssle tshsatnh athnet hCeooCrodoinrdaitnedat eUdniUvnerisvaelr sTailme Tim(UeT(CU-T5C). -O5)f. tOhef ttohtealt ortaailnr gaianugea usgtaetisotnasti oin sthine sthtuedsyt uddoymdaoinm, asiinx ,asriex wariethwini tthhien bthasein, babsient,wbeeetnw 1e7e7n41 m77a4slm (IaNslH(I)N anHd) 3a9n2d03 m92a0slm (SaAslC(S).A DCa)t.aD coatvaercaogve riang tehein sttuhdeys tpuedryiopde rraiondges ranges between 63.8% (INH) and 99.5% (SIC) at rain gauge stations and between 12% (INT) and 36.5% (CHI) at hydrometric stations, except for the PIS station, which has the longest record of hourly streamflows in the basin (98.3%). The distribution of rain gauge and Remote Sens. 2021, 13, 826 4 of 18 hydrometric stations in the study domain is shown in Figure 1. There is a greater density of rain gauge stations towards the northwest of the study domain and a low density towards the southeast. Hydrometric stations are distributed heterogeneously throughout the entire basin, with the middle basin (PIS) having the longest data record. Therefore, of the four stations collected, only one was selected for calibrating the GR4H model, and the remaining three were selected for its validation. Table 1. Detail of main physiographic subbasin characteristics. N◦ Station * Area Mean Elevation Mean Slope Length Subbasin [km2] [masl] [◦] [km] 1 SAL 2042 4753 11 45 2 - 1743 4167 13 41 3 - 290 4317 19 29 4 - 686 4635 18 17 5 - 42 3766 17 9 6 - 1192 4010 17 58 7 PIS 906 3725 17 3 8 - 1113 3858 20 59 9 - 766 3733 13 13 10 CHI 401 4091 23 17 11 INT 411 3791 27 7 (*) Hydrometric stations at the subbasin outlet. Table 2. Selected ground-based stations with hourly records from 1 January 2016 to 31 March 2020. Type Station Abrev. Longitude Latitude Elevation Coverage [ºW] [ºS] [masl] [%] Acjanaco Gore AGR 71.62 13.20 3466.11 85.63 Calca CAL 71.96 13.33 2921.24 94.01 Casaccancha CAS 72.30 13.99 4033.16 84.62 Huayllabamba HUA 72.45 13.27 2976.55 86.86 Intihuatana H INH 72.56 13.17 1774.23 63.86 Intihuatana M INM 72.56 13.17 1778.23 89.81 Machupicchu MAC 72.55 13.18 2399.80 68.18 Pluviometric Marcapata Gore MAR 70.90 13.50 1792.76 83.05 Qorihuayrachina QOR 72.43 13.22 2517.25 96.41 Salcca SAC 71.23 14.17 3920.10 87.75 San Pablo SPB 72.62 13.03 1228.11 90.45 Santa Teresa STR 72.59 13.13 1491.43 64.97 Santo Tomas STM 72.10 14.45 3665.48 95.58 Sicuani SIC 71.24 14.24 3534.95 99.53 Soraypampa SOR 72.57 13.40 3842.32 97.78 Intihuatana km105 INT 72.53 13.18 1774.72 12.01 Hydrometric Chilca CHI 72.34 13.22 2475.28 36.57 Pisac PIS 71.84 13.43 2791.65 98.33 Salcca SAL 71.23 14.17 3918.71 31.39 2.3. Satellite Precipitation Products (SPPs) This study seeks to evaluate the feasibility of implementing a near real-time flood monitoring system at the operational level so that only gridded products with latencies less than 12 h (relative to UTC-5) were selected. Data from the products IMERG-E [30], GSMaP-NRT [31], CMORPH [32], and HE [33] were downloaded and processed for the period from 1 January 2016, to 31 March 2020. The selected SPPs have hourly temporal resolutions and spatial resolutions that vary between 0.1º (~10 km) and 0.05º (~5 km). The details of the characteristics of the selected SPPs are presented in Table 3. Remote Sens. 2021, 13, 826 5 of 18 Table 3. Resume of near real-time satellite precipitation products (SPPs) selected. Product Version Short Name Institution Resolution Latency Integrated Multi-satellite Retrievals for GPM Early V06B IMERG-E NASA 0.1º × 0.1º 5 h Global Satellite Mapping of Precipitation Near Real-Time V06 GSMaP-NRT JAXXA 0.1º × 0.1º 5 h Climate Prediction Center Morphing Method v0.x & v1.0 CMORPH NOAA/CPC 0.08º ×0.08º 8 h HydroEstimator - HE NOAA/NESDIS 0.05º ×0.05º 3 h Note: SPP’s latency is refered to the local time (UTC-5). 3. Methods 3.1. Statistical Evaluation of SPPs against Rain Gauges Because a pixel-to-pixel comparison requires a high density of stations to reduce the uncertainty of the interpolations, added to the different resolutions of the SPPs, in this study, the SPPs were evaluated using a point-to-pixel approach. The pixels of each SPP where the rain gauge stations are located in the study domain were extracted. SPP performance was evaluated through the calculation of statistics (Table 4) that include the relative error (RE), coefficient of correlation (R), root mean square error (RMSE), and mean absolute error (MAE) using the hourly data of all stations for each month of the year. Table 4. Statistical metrics and their corresponding equations used for evaluating the meteorological performance of SPPs. Statistical Metric Unit Equation Optimal Value n Relative Error (RE) - RE ∑i=1(Si−G = i) n 0 ∑i=1 Gi n - ∑i=1[(Si−S)(Gi−G)] 1 Coefficient of Correlation (R) R = √ √ √ n − 2 ∑ n 2 i=1(Si S) ∑i=1(Gi−G) mm/h n 2 0 Root Mean Square Error (RMSE) RMSE = 1 n ∑ (Si − Gi) i=1 n mm/h 1 0 Mean Absolute Error (MAE) MAE = n ∑ |Si − Gi| i=1 Note: n, number of samples; Gi , precipitation from rain gauges; Si , precipitation estimates from satellite products. Likewise, the representation of the diurnal cycle of precipitation in the rainy season (November to March) was evaluated for each station point. For this, only the days with daily precipitation greater than 15 mm were selected, and the coefficient of correlation was calculated for each station point and SPPs. This will allow evaluating not only the magnitude of the SPPs but also their diurnal variability. 3.2. Bias Correction of SPPs In this work, rain gauge observations and the SPPs was merged using a simple bias adjustment correction [34]. For this, in each time interval (t), the pixel corresponding to the station point is extracted, and the bias to the station point (Ybias) is calculated as: Ybiast = Gt − St (1) then Ybias is interpolated by inverse distance weighting (IDW), using a weighting factor of 2, to generate a bias map (Mbias). Finally, the corrected bias is added to the SPPs (SPP′), as shown below: SPP′t = SPPt + Mbiast (2) Remote Sens. 2021, 13, x FOR PEER REVIEW 6 of 18 3.2. Bias correction of SPPs In this work, rain gauge observations and the SPPs was merged using a simple bias adjustment correction [34]. For this, in each time interval (t), the pixel corresponding to the station point is extracted, and the bias to the station point (Ybias) is calculated as: Ybiast = Gt − St (1) then Ybias is interpolated by inverse distance weighting (IDW), using a weighting factor of 2, to generate a bias map (Mbias). Finally, the corrected bias is added to the SPPs (SPP'), as shown below: SPP′t = SPPt + Mbiast (2) 3.3. Semi-distributed GR4H Model 3.3.1. Model Description and Setup TRhemeo teGSenRs.420H21, 1m3, 82o6del [23], which is a variant of the GR4J model of [35] with an6 hofo18urly time step, is used to simulate the runoff in each of the 11 sub-basins of the Vilcanota River (Table 1). The structure of3 t.3h. eSe mGi-RDi4stHrib umtedoGdRe4Hl iMso dsehl own in Figure 2. The contribution of pre- cipitation after interceptio3n.3 .(1P. Mno)d iesl D desicvriipdtieonda nind tSoet udpirect runoff (Pn−Ps) and infiltration (Ps) The GR4H model [23], which is a variant of the GR4J model of [35] with an hourly in the soil moisture reservtiomiers toepf, ims uasexdimto suimul astetothrearugneo ffcianpeaachciotfyth ee1q1usuabl- btaosi nXs o1f.t hTehVielc apnoetracRoivleartion (Perc) from the soil moist(uTarbele 1r)e. sTehrevsotruirct uirse omf tihxeeGdR 4wH imthod etlhise shdoiwrnecint Friugunreo2f.f Tthoe cfonrtrmibu ttihoneo ftotal precipitation after interception (Pn) is divided into direct runoff (Pn−Ps) and infiltration runoff (Pr). Then, the tota(Pl s)riunnthoe sfofi limso pistaurretrietsieorvnoeirdof minaxtiom utmwsoto rragoeucatpinacgity peqruoalcteosXs1.eTsh eupsericnoglat iaon rate defined by parameter X2. (PNericn) ferotym tpheersociel mnot isotufr ethreese rrvuoniroisfmf iixse drwoiuthtethde dbiryec tthruen ocffotmo fobrimntahteitootnal of a runoff (Pr). Then, the total runoff is partitioned into two routing processes using a rate unit hydrograph (UH1) anddefi an endobynplianraemaerte rreX2s.eNrivnoetyirp oerfc emnt aofxtihme runmoff cisarpouatecditby tXhe3c.o Tmhbiena otiothn eofra 10% is routed by another unit hunyidt hryodrgorgarapphh( U(UH1H) an2d).a Tnohnelin beaarsrees eprveoririoofdm aoxifm UumHc2ap a(c2itXy 4X3). iTsh edoothuerb1l0e% that is routed by another unit hydrograph (UH2). The base period of UH2 (2X4) is double that of UH1 (X4). of UH1 (X4). Figure 2. FTighuree s2.tTrhuecsttruurcteu roe fo ftthheeG GR4RH4mHod eml, toakdeenlf,r otmak[3e6]n. from [36]. To implement a flood forecast model in each sub-basin for operational purposes in To implement a floodth efSoErNeAcMasHtI mandogdiveeln ithna tethaecGhR 4sHubm-obdealsisinan faoggrr eogapteemraotdieol,nthaelM pusukrinpgoumse–s in Cunge method [37] was used to transfer the volume of runoff generated in each sub-basin the SENAMHI and given (tsehmait-d tishtreib uGteRd 4mHod eml). oIndtehils iwsa ya,nth eahgogurrleygfloawtew mas osimduella,t etdhfeor M11 ususbk-biansignsum– Cunge method [37] was usaenddr itvoer tsreactniosnfse(Tra tbhlee1 )v. Aodlduitmionea lolyf, trhue nGRo4fHf gmeodneel rreaqtueirdes itnhe einapcuht o sf uhobu-rlbyasin precipitation (PHR) and potential evapotranspiration (ETPHR) data for each sub-basin. (semi-distributed model).T Ihne mtheainsP wHRaoyf ,e atchesu bh-obausirnlwya sflcoalwcu lawtedasfo rsiemachutlimateeindte rfvoarl a1s1th esumbea-nboafsins and river sections (Table 1th)e. bAiasd-cdorirteicotend aSPllPyg,r itdhs,ew hGileRt4heHET mPHoR dwaesle srteimqauteidrebys thtehHea irgnrpeauvets –oSafm haoniurly method [38] from the disaggregation of the daily climatic data (1981–2016) of the mean precipitation (PHR) and potential evapotranspiration (ETPHR) data for each sub-basin. The Remote Sens. 2021, 13, 826 7 of 18 air temperature from the Peruvian Interpolated data of SENAMHI’s Climatological and hydrological Observations (PISCO) product. The PISCO product is a gridded database of precipitation and temperature that coverage throughout the Peruvian territory for the period 1981-2016. The temperature product is generated by applying geostatistical techniques to combining air temperature data from MODIS images and observations from weather stations nationwide and is available online: http://iridl.ldeo.columbia.edu/ SOURCES/.SENAMHI/.HSR/.PISCO/ (accessed on 7 April 2019). Hydrological modeling was performed using the “airGR” package [39] in R language. 3.3.2. Model Calibration, Validation, and Verification The calibration of the model considers the hourly recording at PIS stations from 1 January 2016, to 31 March 2018. The first 300 values were considered as part of the warm- up period to reduce the uncertainty associated with the initial conditions of the model. The parameters of the GR4H model (X1, X2, X3 and X4) were automatically calibrated using the Shuffled Complex Evolutionary algorithm [40], considering the Kling–Gupta efficiency criterion [41] as an objective function. The mean absolute relative error (MARE), percentage bias (PBIAS), and root mean square error (RMSE) were used to evaluate the model performance in terms of low flows representation, model bias, and high flows errors, respectively. The summary of the statistics used to evaluate the performance of the hydrological model is detailed in Table 5. The validation consisted of evaluating the outputs of the GR4H model at the PIS station for the period from 1 April 2018, to 31 March 2020, using the parameters obtained in the calibration. The verification process consisted of evaluating the model performance at the remaining hydrometric stations. For this, the hourly flow records available between 1 August 2018, and 31 March 2020, at stations SAL (upstream of PIS) and CHI and INT (down- stream of PIS) were used. Likewise, this study evaluates the performance of the GR4H model by using the ensemble mean of simulations generated by IMERG-E’, GSMaP-NRT’, CMORPH’, and HE’. In this way, the aim is to evaluate the representativeness of the mean of the simulations commonly used in a flow monitoring system at the operational level. Table 5. Statistical metrics and their corresponding equations used for evaluating the hydrological performance of SPPs. Statistical Metric Unit √ Equation Optimal Value KGE = 1− (r− 1)2 + (α− 1)2 + (β− 1)2 ∑n [(X −X)(O −O)] Kling–Gupta efficiency (KGE) - r = √ i=1 i √ i n − 2 n − 2 1 ∑i=1(Xi X) ∑i=1(Oi O) α = σX µX σ , β = O µO n - 1 |Xi−Oi | 1 Mean Absolute Relative Error (MARE) MARE = 1− n ∑ O i 1 i = n % Percentage Bias (PBIAS) PBIAS =√100 ∑i=1(Xi−Oi) n 0 ∑i=1 Oi m3/s n 2 0 Root Mean Square Error (RMSE) RMSE = 1 n ∑ (Xi −Oi) i=1 Note: n, number of samples; Oi, observed streamflow; Xi, simulated streamflow. 4. Results 4.1. Validation of SPP against Rain Gauges Figure 3 shows the seasonal variation in the metrics calculated to evaluate the meteo- rological performance of the SPPs using all the stations within the study domain. Among all of them, the RE (Figure 3a) shows that the differences between the SPPs are more pro- nounced. Specifically, IMERG-E tends to underestimate hourly precipitation during much of the year, reaching underestimates of up to −53% during the rainy season (November Remote Sens. 2021, 13, 826 8 of 18 Remote Sens. 2021, 13, x FOR PEER REVIEW 8 of 18 to March). In contrast, GSMaP-NRT predominantly overestimates up to 42.5% for the same period. The CMORPH and HE products show moderate differences between the dry 4. Results period (April to October), but they coincide in underestimations close to −20% in the rainy 4.1. Validatipone roifo SdP. PT ahgeaiRnsvt aRlauiens G(aFuiggeusr e 3b) are less than 0.30, as expected for an hourly time scale; Figureh 3o wsheovwesr, tfihtei smeapsronvaels vdauriraintigonth ine wthet mpertrioicds. cIanlctuelramtesdo tfoR eMvaSlEuaatne dthMe mAEet(eF- igure 3c,d), the orological pGeSrfMoramPa-NncReT ofp trhoed SuPcPt sp uressinengt asllt htheeg srteaatitoensts bwiaitshoinf uthpe tsotu1d.1y3 d7oamndain0. 4A3mmomng/ h, respectively, all of them,b tehtew ReEe n(FFigeburue a3ray) sahnodwMs tahracth t.he differences between the SPPs are more pro- nounced. SpecifiAcadlldyi,t IiMonEaRllGy,-EF tiegnudrse to4 usnhdoewresstihmeacteo hrroeularltyio pnre(cRip)itmataiopno dfutrhineg dmiurcnha l cycle of the of the year,o rbesaechrvinegd upnrdeecripesittiamtiaotnesv oefr suups toth −e5S3%PP deusrtiinmg athtees rfaoinryd saeialysopnr e(Nciopviteamtiboenr atob ove 15 mm/d. March). In Tcohnetrsapsat,t iGaSl MdiasPtr-iNbuRtTio pnreodfoRmsinhaonwtlsy goovoedresptoimsiatitvese ufipts to(R 42≥.5%0. 8f)ora tthset astaimone points located period. Thep CreMdoOmRPinHa natnldy HinEt phreosdouucttsh wsheoswt pmoordtieornateo fdtifhfeersetnucdesy bdetowmeaeinn t.heT hderyI MpeE- RG-E product riod (April( tFoi gOucrteob4ear))h, basutt htheehyi gcohiensctidneu minb uenrdoefrpesotiinmtsatwioinths cRlovsae ltuoe −s2g0r%ea itne rththe arnai0n.y8 , including the period. Thes tRa tvioanlusews (iFthiginurteh 3ebb) aasrien l,ewssh tihleanv 0e.r3y0,l oaws evxpaleucteesdo ffoRr a(n0 .h2o