## Abstract

This study aims to develop a stochastic method (SM_GSTR) for generating short-time (i.e., hourly) rainstorm events at all grids (named gridded rainstorm events) in a region. The proposed SM_GSTR model is developed by the non-normal correlated multivariate Monte Carlo simulation (MMCS) method with the statistical properties and spatiotemporal correlation structures of the four event-based gridded rainfall characteristics. The radar-based rainfall data on 20 typhoon events at 336 grids in a basin located in north Taiwan, Nankan River watershed, are used in the model development and demonstration. The results from the model demonstration indicate that the proposed SM_GSTR model can reproduce a great number of gridded rainfall characteristics, of which, the statistical properties in time and space have a good fit to those from the observations in association with the acceptable deviation; thus, it can reasonably emulate the behavior of the rain field in both time and space. It is expected that the resulting massive rainfall-induced disasters (e.g., inundation and landslide) from the physical-based numerical model with the simulated gridded rainstorms by the proposed SM_GSTR model can be applied to establish an alternative artificial intelligence (AI) model for effectively forecasting the hydrologic variables (e.g., runoff and water-level).

## Highlights

A stochastic model for generating a rain field is develop under consideration of correlation of rainfalls in time and space.

The rainstorm event consists of five gridded rainfall characteristics.

The resulting rain field is composed of simulated rainstorm at all grids.

A correlated multivariate Monte Carlo method is applied.

The statistical properties of simulated rainstorms at all grids can be persisted.

## INTRODUCTION

Recently, due to climate change and the occurrence of extreme rainstorm events with high resolution in time and space, investigations related to rainfall-induced disasters have played an essential role in relevant hydrologic applications, such as flood forecasts (e.g., Chang *et al.* 2014; Bhola *et al.* 2018; Zanchetta & Coulibaly 2020), the determination of rainfall thresholds (e.g., Wu *et al.* 2015; Chang *et al.* 2018; Ran *et al.* 2018) and the regulation of the flood risk map (e.g., Amarnath 2014; Doroszkiewicz *et al.* 2019; Aryal *et al.* 2020). For example, the results from numerical simulation show that the inundation depths at the locations and corresponding flooding areas are important and advantageous in relation to the hydrological/hydraulic analyses regarding the quantification and evaluation of the flood-induced damage (Smith 1997). However, it remains generally difficult to measure the real-time practical inundation-depth hydrograph and even the induced flooding area due to the limitation of measurement equipment or hindrance in data acquisition, processing, and analysis (e.g., Amarnath 2014; Park *et al.* 2018); this situation possibly leads to a lack of reliable references regarding the fine spatiotemporal resolution required for the parameter calibration of the hydraulic numerical models (Shen *et al.* 2019). Despite the fact that the possible maximum flooding area caused by a rainstorm event can be estimated in accordance with the observed maximum inundation-depth recorded at particular locations used to calibrate the parameters of a hydraulic numerical simulation model, the resulting temporal and spatial hydrological features (e.g., inundation depth) include uncertainties with high likelihood attributed to insufficient long-term observations (Hardesty *et al.* 2018). In order to successfully carry out high-resolution inundation simulation in time and space, the rainfall feature with high spatiotemporal resolution (i.e., gridded precipitation or high-density gauged rainfall) should be taken into account (e.g., Lovat *et al.* 2019; Try *et al.* 2020; Zanchetta & Coulibaly 2020). Although improving the spatial density of the rain-gauge networks can effectively capture the rainfall features representing the components of the hydrological cycle, the high cost of equipment and associated maintenance require consideration (Terink *et al.* 2018). Therefore, rain-field simulation modeling is required in implementing flood and inundation simulations through the hydrologic/hydraulic model (e.g., Wheater *et al.* 2005; Paschalis *et al.* 2013; Jiang & Tatano 2015; Wu *et al.* 2015; Chang *et al.* 2018; Zanchetta & Coulibaly 2020).

Numerous methods for reproducing rain fields are developed on the basis of physically based models and stochastic modeling, respectively; the commonly used physically based models are the general circulation models (GCMs) which predominantly simulate and predict climatic features (e.g., precipitation and temperature) on a global scale (e.g., Reed 1986; Moron *et al.* 2008; Khan *et al.* 2018; Ahmed *et al.* 2019; Liu *et al.* 2020). Nevertheless, GCMs commonly generate gridded precipitation based on the project-climate conditions with large uncertainties caused by the concept assumption, model structure, and calibrated parameters (e.g., Yip *et al.* 2011; Khan *et al.* 2018).

Regarding the stochastic approaches used in rainfall simulation, the statistical moments of the point-based precipitation (i.e., gridded rainfall characteristics) are necessary to reproduce rainfall processes (e.g., Wheater *et al.* 2005; Yang *et al.* 2005; Li *et al.* 2012; Leblois & Crutin 2013; Paschalis *et al.* 2013; Jiang & Tatano 2015; Gao *et al.* 2018; Park *et al.* 2018; Qiang *et al.* 2020), For example, Wheater *et al.* (2005) utilized the generalized linear models (GLMs) to generate the daily rainfall using the corresponding spatial structure extracted from observations; they are then combined with the estimated hourly based sequence by a single-site Poisson process under the assumption of temporal stationarity, enabling the simulated hourly rainfall series at different locations in a region to be obtained. Leblois & Crutin (2013) presented an approach for the simulation of the rainfall-field sequences within the advection field through a classical Gaussian regional simulation technique and the turning-band method. Paschalis *et al.* (2013) developed a stochastic model for simulating high-resolution spatiotemporal rainfall (named STREAP) by means of the multivariate stochastic simulation method with the spatial and temporal statistical structure of observed radar-derived rainfall data under the assumption of isotropicity in space (i.e., the identical rainfall accumulation for the entire spatial domain) and the isotropic-simulated rainfall fields. Yang *et al.* (2005) adopted the observations of gauged daily rainfall to identify the best-fit GLM and to quantify the statistical properties, including the spatial correlation structure regarding daily rainfall, the joint distribution for the wet-day pattern of rainfall occurrence and the probability distribution of daily rainfall amount at wet sites used in the rain-field simulation. Gao *et al.* (2018) presented a daily rainfall simulation model by producing the rainfall events comprising simulated rainfall characteristics in terms of rainfall duration, storm depth, and dimensionless storm pattern defined by Wu *et al.* (2006) at gauged sites. In Gao *et al.*’s method, the generation of daily rainfall is implemented through the non-normal correlated multivariate Monte Carlo simulation (MMCS) approach with their copula dependence functions, developed based on the Nataf transformation in association with a linear correlation matrix of random variables (Lebrun & Dutrfoy 2009). By hypothesizing isotoropicity in rainfall characteristics within a minor region, the resulting rainfall events exhibit similar behavior in time to those at the gauges. Qiang *et al.* (2020) proposed a rainstorm generator for reproducing the rainfall process at various locations by calculating the areal rain field using the interpolation approach with gauged simulated rainfall under consideration of the stochastic dependence matrix in space in terms of a copula function (Jiang & Tatano 2015; Gao *et al.* 2018); furthermore, Qiang's rainstorm generator can be extended to simulate rain fields in response to climate change by adjusting the parameters of the marginal probability distributions. In addition to stochastic modeling with the Monte Carlo simulation, Park *et al.* (2018) proposed a hybrid stochastic rainfall model to generate rainfall characteristics. In Park *et al.*’s model, the monthly rainfalls are first simulated using the seasonal autoregressive integrated moving average approach, and the monthly rainfalls are disaggregated into the hourly based rainfalls through the modified Bartlett-Lewis rectangular pulse (MBLRP) model. Apart from the simulation models with prior statistical properties of rainfall features, Oriani *et al.* (2014) developed a nonparametric approach to generate the rainfall sequence through direct sampling technique without taking into account the aforementioned prior statistics. Origni's method can provide the time series of daily rainfall and the associated temporal patterns by means of the autoregressive and Markov chain-based methods on the basis of the training image obtained from the previous rainfall record; however, it has inherent limitations in representing extremes (Li *et al.* 2012).

Of the numerous rainfall-simulation models developed, most focus on simulating the correlated rainfall amount/depth for a specific rainfall duration at various locations in the region of interest under an assumption of the isotropic storm pattern in terms of the identical rainfall distribution in time (e.g., Paschalis *et al.* 2013; Jiang & Tatano 2015; Yang *et al.* 2005; Gao *et al.* 2018; Qiang *et al.* 2020). That is to say, the past investigations mainly emphasized the assessment regarding the effects of temporal and spatial changes in rainfall on the relevant hydrologic and hydraulic analysis. However, in reality, in addition to the spatial or temporal uncertainties in the rainfall, the spatiotemporal resolution and variation in the rain field possibly influences the urban/watershed hydrologic responses required in flood-hydrological modeling and the relevant induced disaster prevention and mitigation (e.g., Reed *et al.* 2007; Leblois & Crutin 2013; Jiang & Tatano 2015; Park *et al.* 2018). Therefore, this study aims to develop a stochastic model for generating gridded short-term (i.e., hourly) rainstorm events (called SM_GSTR model) by taking into account the spatiotemporal dependence in the rain field. It is expected that hydraulic dynamic modeling with the resulting rainstorms at all grids (i.e., gridded rainstorms) from the SM_GSTR model can provide more detailed information on rainfall-induced forecasts with fully high resolution in time and space for the decision-making and infrastructure-planning of disaster early warning and proofing systems.

## METHODOLOGY

### Model concept

This study intends to develop a stochastic model for generating a large number of gridded rainstorm events, named the SM_GSTR model. In detail, the proposed SM_GSTR model is established based on Wu's method (Wu *et al.* 2006) involving a stochastic model for generating event-based rainstorm events comprising the at-site rainfall characteristics, including rainfall duration, storm depth, and storm pattern without considering their spatial correlations among rainfall characteristics. Also, the spatial heterogeneity regarding rainfall significantly impacts the shape of the hydrological response, especially for the combination of catchment features and rainfall characteristics (Volpi *et al.* 2012), By doing so, this study modifies Wu's model by taking into account both spatial and temporal correlations regarding rainfall characteristics extracted from the gridded rainfall data (called gridded rainfall characteristics). It is noticed that the main difference of the proposed SM_GSTR model from Wu's method is the correlation structure in both time and space considered in the generation of rainstorm events at all grids within a region.

In developing the proposed SM_GSTR model in accordance with the model framework, the observations of the rainfall characteristics at all grids within the study area, defined as gridded rainfall characteristics, are extracted from historical hourly rainfall data. Uncertainty analysis is then carried out to quantify the corresponding statistical properties regarding time and space, including the first four statistical moments and correlations; in so doing, the resulting statistical properties of gridded rainfall characteristics are used to reproduce the short-term rainstorm events at all grids in a region, composed of the simulated gridded rainfall characteristics by the non-normal correlated MMCS method (Wu *et al.* 2006). The aforementioned relevant methods and concepts are addressed below.

### Brief introduction to stochastic generation of at-site rainstorms

#### At-site rainfall characteristics

*et al.*(2006), in which the hourly rainstorm events at a single rain gauge can be characterized into three characteristics in advance, the inter-event time, rainfall duration, rainfall amount and storm pattern (see Figure 1). Among the at-site rainfall characteristics, the storm pattern comprises the dimensionless rainfall at various dimensionless time steps obtained by the following equation (Wu

*et al.*2006):in which ∈[0, 1] represents the dimensionless time;

*d*is the storm duration; ∈[0, 1] stands for the cumulative dimensionless rainfall; is the incremental dimensionless rainfall amount;

*D*is the cumulative rainfall depth at time

_{t}*t*; and

*D*denotes the total rainfall depth.

_{d}#### Simulation of at-site rainfall duration and storm depth

*et al.*(1997) is adopted herein to simultaneously simulate the rainfall duration and storm depth for a specific number of simulated events given in advance. The non-normal correlated MMCS method has the ability to preserve the marginal distributions of random variables and their correlation structure without requiring the complete joint distribution through three steps: transformation to standard normal space, orthogonal transform, and inverse transformation, which can be implemented by the Nataf bivariate distribution model (Nataf 1962):where

*X*and

_{i}*X*are the correlated variables at the points

_{j}*i*and

*j*, respectively, with the means and , the standard deviations and , and the correlation coefficient ;

*Z*and

_{i}*Z*are corresponding bivariate standard normal variables to the variable

_{j}*X*and

_{i}*X*with the correlation coefficient and the joint standard normal density function

_{j}*ϕ*(·). In referring to Equation (2), the MMCS method can generate multivariate with high correlation in time or space on the basis of the selected marginal distributions and correlation coefficients (Chang

_{ij}*et al.*1997).

*T*, depending on the marginal distributions and correlation of

_{ij}*X*and

_{i}*X*, can be determined to modify in the original space to in the normal space by:

_{j}*et al.*(1997) derived a number of formulae for calculating the transformation factor

*T*between the correlated multivariate standard normal variables and uncorrelated standard normal variables by which random variates can be easily generated. Eventually, upon the generation of multivariate standard normal random variates , one can obtain the corresponding random variates in the original space by the following inverse transformation:where is the marginal cumulative distribution function (CDF) of random variable X

_{ij}_{i}in the original space, and Ф(·) is the standard normal CDF. The framework of generating the correlated multivariate by means of the MMCS method can be referred to in Figure 2.

#### Simulation of storm pattern

In Wu's method, the rainfall duration and storm depth regarded as the spatial rainfall characteristics can be directly generated by the MMCS method. However, the storm pattern is defined as the time distribution of rainfall, comprising dimensionless rainfall in which the time-axis is divided into 12 increments (Wu *et al.* 2006). Accordingly, two special constraints of dimensionless rainstorm patterns must be observed: (1) the summation of dimensionless rainfall must be equal to 1: P_{1/12} + P_{2/12} + …+ P_{12/12} = 1 as 0 ≤ F_{1/12} ≤ F_{2/12} ≤ …≤ F_{12/12} = 1; and (2) the dimensionless rainfall must be between 0 and 1: 0 ≤ P*τ* ≤ 1 for *τ* =1/12, 2/12, …, 12/12.

*et al.*(2006) utilized the Johnson distribution system, which is more flexible in fitting the rainstorm pattern than are other distributions (Fang & Tung 1996), to standardize the log-ratio variates instead of the standard normalization equation, i.e., Equation (6), used in the non-normal correlated MMCS with Equation (2). It is well-known that the Johnson distribution is a system of frequency curve consisting of four parameters (Johnson 1949):where Z is standard normal random variable; g(.) is a monotonic function of the original random variable

*X*; and are the location and scale parameters, respectively. In general, the Johnson distribution system is grouped into three types:

Of the three types of Johnson distribution, there is very good agreement between the measured and estimated rain variables using the Johnson distribution system, especially for the bounded system () (D'Adderio *et al.* 2016). Hill *et al.* (1976) presented an algorithm to identify the best-fit type of the Johnson distribution system and calibrate the associated parameters listed in Equation (10) by means of matching the first four product-moments of the random variables.

Consequently, the framework of generating at-site rainfall characteristics through Wu's method includes five steps: (1) quantification of statistical properties of at-site rainfall characteristics; (2) generation of rainfall duration and storm depth; (3) simulation of storm patterns composed of generated dimensionless rainfalls; and (4) reproduction of at-site rainstorms by combining the simulated at-site rainfall characteristics as shown in Figure 3.

### Development of the proposed SM_GSTR model

#### Definition of gridded rainfall characteristics

*et al.*2015; Cristiano

*et al.*2017; Zhu

*et al.*2018; Qiang

*et al.*2020). Hence, to consider the spatiotemporal correlation among at-site rainfall characteristics, the SM_GSTR model treats the rainfall characteristics at various grids (i.e., gridded rainfall characteristics) within the watershed as the spatially and temporally correlated variables. In detail, the gridded rainfall depths can be treated as the spatially correlated variables, whereas the storm patterns in terms of dimensionless rainfall at entire grids should not only be the temporal variable, but also the spatial one. Specifically, according to hydrologic isotropicity, the rainfalls at different locations during a rainstorm event are supposed to exhibit similar temporal changes in rainfall, especially for the minor basin (Cristiano

*et al.*2017); however, in the case of the high variability of rainfall in space and time, the rainfall pattern in time and space should be associated with somewhat differences/bias. By doing so, in this study, in response to the effect of uncertainty of rainfall in time and space, the storm pattern at all grids (i.e., gridded storm pattern) within a watershed could be separated into two components: the areal average dimensionless rainfall and the associated difference from the areal average (called the gridded deviation) could be calculated by the following equations:in which and stand for the dimensionless rainfall at the dimensionless time k

^{th}and resulting areal average for the j

^{th}grid; is the cumulative dimensionless rainfall at the dimensionless time for the j

^{th}grid; is the areal average of the dimensionless rainfall at the dimensionless ; and serves as the difference in the cumulative dimensionless rainfall at the dimensionless time for the j

^{th}grid from the corresponding areal average. Note that the areal average is calculated through the following equation:where serves as the number of grids in the steady area; stands for the cumulative dimensionless rainfall for the dimensionless time at the j

^{th}grid. Since the and can be treated as the temporal and spatiotemporal variates, respectively, the deviation is accordingly regarded as the spatiotemporal variable.

In summary, this study defines gridded short-term rainstorms in terms of gridded rainfall characteristics, including the rainfall duration of events, the rainfall depths at all grids (i.e., gridded rainfall depths), and the areal average storm pattern consisting of cumulative dimensionless rainfalls and the associated deviation at different dimensionless times. Figure 4 presents the process of extracting the aforementioned gridded rainfall characteristics from the observed gridded hyetographs regarding historical rainstorms.

#### Stochastic generation of gridded rainfall characteristics

With respect to Wu's method, the non-normal correlated MMCS method (Wu *et al.* 2006) is utilized to generate the highly temporally correlated rainfall characteristics. By doing so, the proposed SM_GSTR model also applies the MMCS method to simulate a number of the gridded rainfall characteristics in response to the uncertainties in time and space. In detail, the process of generating the gridded rainfall characteristics by means of the proposed SM_GSTR model can be grouped into four parts: the simulation of the rainfall duration and rainfall depths at all grids (i.e., gridded rainfall depths), the simulation of the areal average storm patterns, and the generation of the deviation of gridded storm pattern. The detailed concept of the aforementioned generation process is introduced in the following:

- 1)
Rainfall duration and gridded rainfall depths: Regarding the proposed SM_GSTR model, during a simulated rainstorm event, the corresponding rainfall duration and a number of the rainfall depths at the grids of interest (i.e., gridded rainfall depths) are correlated variables, which differ from Wu's method for simulating at-site events. Therefore, in the case of N grids used in the field, the MMCS model is adopted herein to simultaneously simulate (

*N*+ 1) variables, including rainfall duration and storm depths of the N grids. - 2)

*τ*for the j

^{th}grid; denotes the simulated areal average of cumulative dimensionless rainfall at the dimensionless time

*τ*and represents the generated deviation of the cumulative dimensionless rainfall at the dimensionless time

*τ*for the j

^{th}grid. However, in the case of the unreasonable simulated cumulative dimensionless rainfall of storm pattern obtained through Equation (15), i.e., , this study modifies the through the optimization method with the flooding equation to find out an adjusting factor of the gridded mean deviation:in which and stand for the modified simulation of cumulative dimensionless and the corresponding adjusting factor, respectively. In Equation (16), the adjusting factor of the gridded deviation () can be determined by the optimization method, such as the genetic algorithm (GA), subject to a condition as:

Eventually, using a great number of simulations of the five gridded rainfall characteristics, the resulting hyetographs of simulated rainstorm events at all grids in a region of interest can be obtained through the graphical process as shown in Figure 5.

### Framework of modeling stochastic generation of gridded rainstorms

To sum up the concepts introduced above, the proposed SM_GSTR model includes (1) the characterization of gridded short-term rainstorm events; (2) the quantification of uncertainties of gridded rainstorms in time and space; and (3) the stochastic generation of spatiotemporally correlated rainfall characteristics. Note that in generating the gridded rainfall characteristics, especially for the rainfall durations and gridded rainfall depths, via the proposed SM_GSTR model, the conditions of non-negative gridded rainfall characteristics should be given in order to avoid the unreasonable simulations, i.e., the negative rainfall durations, rainfall depths, and dimensionless rainfalls of storm patterns. The detailed framework of developing the proposed SM_GSTR model can be addressed (see Figure 6) as follows:

Step [1]: Collect the grid-based precipitation of numerous rainstorm events in the study area and extract their gridded rainfall characteristics, i.e., rainfall duration, gridded rainfall depth, areal average of cumulative dimensionless rainfall and the associated deviation at various dimensionless times.

Step [2]: Implement uncertainty analysis for the gridded rainfall characteristics to quantify the uncertainties in terms of the statistical moments and correlation coefficients in time and space.

Step [3]: Generate the event-based rainfall duration and the gridded rainfall depths using the non-normal correlated MMCS method based on their statistical properties and correlation coefficients in time and space.

Step [4]: Generate the areal average of the dimensionless rainfalls based on the statistical moments and their temporal correlation coefficients; then, cumulate the simulated areal averages as the areal average cumulative dimensionless rainfalls.

Step [5]: Simulate the deviation of cumulative dimensionless rainfall at all grids (i.e., gridded deviation) for each dimensionless time based on the corresponding statistical moments and correlation coefficients.

Step [6]: Combine the simulated average of the cumulative dimensionless rainfall at each dimensionless time at Step [4] with the corresponding generated grid-based deviation obtained at Step [5] as the cumulative dimensionless rainfall at each dimensionless time at all grids in the study area.

Step [7]: Combine the resulting simulations of gridded rainfall characteristics from the above steps as the gridded short-term rainstorm events. The graphical process of establishing the proposed SM_GSTR model can be seen in Figure 6.

## STUDY AREA AND DATA

### Introduction to study area

The study area, the Nankan River, in which the corresponding length and the drainage area are 31 km and 224 km^{2}, respectively, is in Taoyuan County of northern Taiwan (see Figure 7); also, on average, the slope approximates 0.0077. In addition, it flows through the Guishan, Taoyuan, and Luzhou Districts in Taoyuan City. Within the Nankan River watershed, there are six riverside parks and three branches, Dongmon Creek, KengZi Creek, and Qiedong Creek. Of these branches, the flooding risk at the Dongmon Creek, previously inundated, has been noticeably reduced by increasing the number of the detention ponds and the cross-section area in the river channel. Additionally, it is well-known that Taiwan Central Weather Bureau (CWB) can provide the quantitative precipitation estimation (QPE) with the spatial resolution of 1.5 km × 1.5 km; hence, there are 336 grids (called QPE grids) located in Nankan River watershed as shown in Figure 7 and the corresponding gridded hourly rainfall data can be utilized in the development and demonstration of the proposed SM_GSTR model.

### Spatiotemporal analysis of the study data

Since the purpose of this study is to develop a stochastic model for simulating a large number of rainstorms with high resolution in time and space, the hourly rainfall datal of 20 rainstorm events at 336 grids within the study area (Nankan River watershed) (see Table 1) are selected as the study data.

Event . | Beginning time . | Ending time . | Duration (hours) . | Areal average rainfall depth (mm) . |
---|---|---|---|---|

1 | 2005/7/17 10:00 | 2005/7/18 18:00 | 33 | 139.8 |

2 | 2005/8/4 09:00 | 2005/8/6 00:00 | 40 | 230.6 |

3 | 2005/8/31 07:00 | 2005/9/1 08:00 | 26 | 125.5 |

4 | 2005/10/1 23:00 | 2005/10/2 18:00 | 20 | 54.5 |

5 | 2008/7/26 22:00 | 2008/7/28 22:00 | 49 | 105.8 |

6 | 2008/9/12 02:00 | 2008/9/15 00:00 | 71 | 403.2 |

7 | 2010/10/21 00:00 | 2010/10/22 13:00 | 38 | 120.8 |

8 | 2012/6/11 21:00 | 2012/6/12 23:00 | 27 | 366.7 |

9 | 2012/6/14 11:00 | 2012/6/15 06:00 | 20 | 18.5 |

10 | 2012/8/26 09:00 | 2012/8/27 08:00 | 24 | 29.8 |

11 | 2013/5/11 01:00 | 2013/5/13 01:00 | 49 | 184.8 |

12 | 2013/7/12 16:00 | 2013/7/13 14:00 | 23 | 75.4 |

13 | 2013/10/4 14:00 | 2013/10/6 23:00 | 58 | 56.3 |

14 | 2014/5/20 20:00 | 2014/5/22 00:00 | 29 | 195.9 |

15 | 2014/7/22 21:00 | 2014/7/24 03:00 | 31 | 43.4 |

16 | 2014/9/21 16:00 | 2014/9/22 12:00 | 21 | 83.5 |

17 | 2015/8/7 08:00 | 2015/8/8 13:00 | 30 | 184.9 |

18 | 2015/9/27 14:00 | 2015/9/29 05:00 | 40 | 138.3 |

19 | 2016/9/26 10:00 | 2016/9/28 04:00 | 43 | 106.3 |

20 | 2017/6/2 10:00 | 2017/6/4 04:00 | 39 | 271.2 |

Event . | Beginning time . | Ending time . | Duration (hours) . | Areal average rainfall depth (mm) . |
---|---|---|---|---|

1 | 2005/7/17 10:00 | 2005/7/18 18:00 | 33 | 139.8 |

2 | 2005/8/4 09:00 | 2005/8/6 00:00 | 40 | 230.6 |

3 | 2005/8/31 07:00 | 2005/9/1 08:00 | 26 | 125.5 |

4 | 2005/10/1 23:00 | 2005/10/2 18:00 | 20 | 54.5 |

5 | 2008/7/26 22:00 | 2008/7/28 22:00 | 49 | 105.8 |

6 | 2008/9/12 02:00 | 2008/9/15 00:00 | 71 | 403.2 |

7 | 2010/10/21 00:00 | 2010/10/22 13:00 | 38 | 120.8 |

8 | 2012/6/11 21:00 | 2012/6/12 23:00 | 27 | 366.7 |

9 | 2012/6/14 11:00 | 2012/6/15 06:00 | 20 | 18.5 |

10 | 2012/8/26 09:00 | 2012/8/27 08:00 | 24 | 29.8 |

11 | 2013/5/11 01:00 | 2013/5/13 01:00 | 49 | 184.8 |

12 | 2013/7/12 16:00 | 2013/7/13 14:00 | 23 | 75.4 |

13 | 2013/10/4 14:00 | 2013/10/6 23:00 | 58 | 56.3 |

14 | 2014/5/20 20:00 | 2014/5/22 00:00 | 29 | 195.9 |

15 | 2014/7/22 21:00 | 2014/7/24 03:00 | 31 | 43.4 |

16 | 2014/9/21 16:00 | 2014/9/22 12:00 | 21 | 83.5 |

17 | 2015/8/7 08:00 | 2015/8/8 13:00 | 30 | 184.9 |

18 | 2015/9/27 14:00 | 2015/9/29 05:00 | 40 | 138.3 |

19 | 2016/9/26 10:00 | 2016/9/28 04:00 | 43 | 106.3 |

20 | 2017/6/2 10:00 | 2017/6/4 04:00 | 39 | 271.2 |

#### Event-based rainfall duration and storm depth

Figure 8 represents the hyetographs of the 20 selected radar-based rainstorm events (2005–2017) provided by Taiwan Central Weather Bureau. In view of Table 1 and Figure 7, the rainfall duration ranges from 20 hours to 70 hours; and, on average, the gridded rainfall depths are located between 40 and 410 mm, implying that the rainstorm events with various rainfall intensities can be used to evaluate the effect of rainstorm scale on the hydrologic response. Also, Figure 9 shows the gridded rainfall depths of 20 rainstorms, indicating that the spatial distribution of the gridded rainfall depths exhibits an obvious change with the rainstorm event considered. For example, in regard to most of the rainstorm events (e.g., EV3, EV13, and EV18), their rainfall depths at all grids gradually approach two centers; namely, the maximum rainfall depths of most rainstorm events can be found approximately near two grids; it can also be found that the rainfall depth for rainstorms (e.g., EV10, EV13, and EV15) at all grids are distributed normally in the study area. The above results indicate significant spatial uncertainty and variation in the gridded rainfall depth in the study area.

#### Gridded storm patterns and associated deviation

In this study, the gridded storm pattern consisting of the cumulative dimensionless rainfalls at various dimensionless times is grouped into two components: the areal average storm pattern and associated deviation at all grids (i.e., gridded deviation) defined as the gridded rainfall characteristics. According to the framework of developing the proposed SM_GSTR model, the uncertainty quantification should be carried out to quantify the corresponding statistical properties in advance (see Figures 10 and 11). From Figure 10, it can be seen that there are obviously different varieties of the areal average storm patterns in the study area; hence, the temporal variations in the storm pattern should be considered. Also, Figure 11 presents the spatial varying trend of gridded deviation regarding the storm pattern. It indicates that the gridded deviation exhibits an obvious varying trend in space; namely, its spatial distribution changes with the dimensionless times and rainstorms. Accordingly, the gridded deviation of storm pattern possibly involves uncertainty in time and space.

To sum up the above results, the gridded rainfall characteristics could be regarded as highly spatially and temporally correlated variables in association with the uncertainties of different degrees. By doing so, as simulating the rain field with high resolution in time and space, it is essential to take into account the uncertainties of gridded rainfall characteristics and their dependence structure in time and space.

## RESULTS AND DISCUSSION

To consider the spatial and temporal variations and correlations in the rain field, a stochastic modeling for generating a great number of gridded rainstorms is established based on the model development framework (see Figure 6). In detail, the uncertainties in the gridded rainfall characteristics extracted from the radar-based rainfall data (i.e., QPE) in the section ‘Introduction to study area’ can be quantified to obtain the spatial and temporal statistical properties; massive gridded rainfall characteristics can then be reproduced through the MMCS method. Eventually, a great number of gridded rainstorms involved in the rain field can be obtained by integrating the gridded rainfall characteristics. The relevant results are evaluated and discussed below.

### Uncertainty analysis for grid-based rainfall characteristics

#### Event-based rainfall and gridded rainfall depth

Figure 12 shows the statistical properties of the rainfall duration and gridded rainfall depth; on average, the durations of gridded rainstorm events frequently occurring in the Nankan River watershed is about 35 hours; whereas the lower and upper bounds of 95% confidence interval are 20 hours and 53 hours, respectively, indicating that 1-day and 2-day events frequently take place in the study area, the Nankan River watershed. In the case of the statistical properties of gridded rainfall depth, its average and standard deviation ranges between 80 and 220 mm as well as 82 and 177 mm. Although its lower bound of 95% confidence interval approximates 20 mm, its upper bound has a wide range from 284 to 336 mm, respectively. This implies that even for the same rainstorm event, the variation in rainfall depth at all grids in the study area obviously varies with the locations of the grids.

#### Areal average of storm pattern and gridded deviation

For the statistical properties of cumulative dimensionless rainfall of the areal average storm pattern as shown in Figure 13 and Table 2, a wider 95% confidence interval can represent that various varieties of storm patterns (e.g., central, advanced, and delayed-shaped patterns) take place in the study area. Also, Table 2 lists the correlation coefficients among the dimensionless rainfall (about −0.3 to 0.97); it shows that the areal averaged cumulative dimensionless rainfalls are strongly correlated with different dimensionless times. In addition, the corresponding statistics of its deviation at each grid should be considered, as shown in Figure 13. Moreover, Figure 3 illustrates the mean and standard deviation of the gridded deviation of storm patterns for the dimensionless time selected (, , , and ). It can be known that the aforementioned statistical moments significantly change with the dimensionless time in which the large moments can be found at the central dimensionless rainfall (_{5}); whereas, the cumulative dimensionless rainfall at the beginning (_{3}) and ending time steps (_{9}) are obviously smaller than ones at the remaining time steps. Also, since the correlation coefficients among the gridded deviations of cumulative dimensionless rainfalls range approximately from −0.87 to 0.95, the gridded deviation should be treated as the spatially correlated variables. By doing so, in the proposed SM_GSTR model, the areal average and gridded deviation of storm patterns must be generated under consideration of the correlation coefficients in time and space given in advance.

Dimensionless time τ . | Pτ_{1}
. | Pτ_{2}
. | Pτ_{3}
. | Pτ_{4}
. | Pτ_{5}
. | Pτ_{6}
. | Pτ_{7}
. | Pτ_{8}
. | Pτ_{9}
. | Pτ_{10}
. | Pτ_{11}
. | Pτ_{12}
. |
---|---|---|---|---|---|---|---|---|---|---|---|---|

Pτ_{1} | 1 | |||||||||||

Pτ_{2} | 0.016 | 1 | ||||||||||

Pτ_{3} | −0.2 | 0.87 | 1 | |||||||||

Pτ_{4} | −0.243 | 0.792 | 0.964 | 1 | ||||||||

Pτ_{5} | −0.244 | 0.749 | 0.93 | 0.966 | 1 | |||||||

Pτ_{6} | −0.231 | 0.736 | 0.907 | 0.938 | 0.964 | 1 | ||||||

Pτ_{7} | −0.187 | 0.646 | 0.804 | 0.831 | 0.847 | 0.88 | 1 | |||||

Pτ_{8} | −0.205 | 0.481 | 0.642 | 0.687 | 0.702 | 0.722 | 0.807 | 1 | ||||

Pτ_{9} | −0.142 | 0.322 | 0.417 | 0.451 | 0.459 | 0.474 | 0.553 | 0.702 | 1 | |||

Pτ_{10} | −0.209 | 0.186 | 0.282 | 0.288 | 0.273 | 0.28 | 0.275 | 0.269 | 0.458 | 1 | ||

Pτ_{11} | −0.176 | −0.037 | 0.046 | 0.028 | 0.003 | −0.029 | −0.103 | −0.211 | −0.102 | 0.076 | 1 | |

Pτ_{12} | 0.02 | 0.017 | 0.015 | 0.013 | 0.021 | 0.026 | 0.035 | 0.008 | −0.015 | 0.007 | −0.006 | 1 |

Dimensionless time τ . | Pτ_{1}
. | Pτ_{2}
. | Pτ_{3}
. | Pτ_{4}
. | Pτ_{5}
. | Pτ_{6}
. | Pτ_{7}
. | Pτ_{8}
. | Pτ_{9}
. | Pτ_{10}
. | Pτ_{11}
. | Pτ_{12}
. |
---|---|---|---|---|---|---|---|---|---|---|---|---|

Pτ_{1} | 1 | |||||||||||

Pτ_{2} | 0.016 | 1 | ||||||||||

Pτ_{3} | −0.2 | 0.87 | 1 | |||||||||

Pτ_{4} | −0.243 | 0.792 | 0.964 | 1 | ||||||||

Pτ_{5} | −0.244 | 0.749 | 0.93 | 0.966 | 1 | |||||||

Pτ_{6} | −0.231 | 0.736 | 0.907 | 0.938 | 0.964 | 1 | ||||||

Pτ_{7} | −0.187 | 0.646 | 0.804 | 0.831 | 0.847 | 0.88 | 1 | |||||

Pτ_{8} | −0.205 | 0.481 | 0.642 | 0.687 | 0.702 | 0.722 | 0.807 | 1 | ||||

Pτ_{9} | −0.142 | 0.322 | 0.417 | 0.451 | 0.459 | 0.474 | 0.553 | 0.702 | 1 | |||

Pτ_{10} | −0.209 | 0.186 | 0.282 | 0.288 | 0.273 | 0.28 | 0.275 | 0.269 | 0.458 | 1 | ||

Pτ_{11} | −0.176 | −0.037 | 0.046 | 0.028 | 0.003 | −0.029 | −0.103 | −0.211 | −0.102 | 0.076 | 1 | |

Pτ_{12} | 0.02 | 0.017 | 0.015 | 0.013 | 0.021 | 0.026 | 0.035 | 0.008 | −0.015 | 0.007 | −0.006 | 1 |

### Evaluation of simulations of gridded rainfall characteristics

Using the proposed SM_GSTR method with the spatial and temporal statistical properties of gridded rainfall characteristics, the resulting 1,000 simulations can be obtained as shown in Figures 14 and 15. To evaluate the applicability of the proposed SM_GSTR method for simulating the gridded rainstorms, the statistical features calculated from 1,000 simulation at all grids within the study area are compared with the results from the observations (see Figures 16–21), including the mean, standard deviation, and 95% confidence interval. In theory, since the gridded rainfall characteristics in time and space are generated by using the non-normal correlated MMCS based on the Nataf-based algorithm, in which the correlation structure among simulations has been proven to be consistent with those from observations by Chang *et al.* (1997), the consistency of the correlation structures in time are evident regarding the gridded rainfall characteristics.

In this study, since the gridded rainfall events are reproduced by combining the simulated gridded rainfall characteristics, this study focuses on evaluating the performance in comparison to the correlation with the observed data by calculating the correlation coefficients among the statistical properties, including the mean, standard deviation, and 95% confidence interval (see Figure 22). In detail, the spatial correlation coefficients (i.e., the gridded rainfall depths and the gridded deviation of the cumulative dimensionless rainfall) and temporal correlation coefficients (i.e., the areal average of cumulative dimensionless rainfalls) are adopted to assess the consistency of correlations regarding the gridded rainfall characteristics. The relevant results and discussion are addressed below.

#### Rainfall duration and gridded rainfall depths

By observing Figure 16, the rainfall durations of 1,000 simulated gridded rainstorm events change from 3 hours to 81 hours. As compared with the statistics from observations, the average of the simulated duration (36 hours) approximates that from the observation (37 hours); and the standard deviation of the simulated duration (12 hours) is nearly equal to that from the observations (13 hours). It can be seen that the lower bound of the 95% confidence interval of the simulated duration (12 hours) is smaller than that from observation (20 hours); accordingly, a larger upper bound of 65 hours, than that from observation (53 hours), should be obtained. This difference implies that there noticeably exists variation in event-based rainfall duration; thus, the rainstorm events of different time scales are applied in the development of the proposed SM_GSTR model in response to the short-term effect on the simulation of gridded rainstorms.

In referring to Figure 17, it can be seen that the mean values of simulated gridded rainfall depth (on average 141.3 mm) is close to that from observation (about 164.2 mm); whereas the standard deviation of simulated gridded rainfall depth (41.2 mm–160.8 mm) has an acceptable relative difference (−0.47) to observed values (78.8 mm–176.8 mm). The lower bound of the simulated gridded rainfall depth (on average 74.1 mm) significantly departs from the observed values (15.7 mm), exhibiting good agreement between the simulations and observations in space; on the contrary, the upper bound (234.2 mm–781.4 mm) is slightly greater than the observations (267.7 mm–686.6 mm). Also, similar to the above statistical features, the change in the upper bound of simulated gridded rainfall depth approximately matches those from the observation.

As a result of the gridded rainfall depths being the spatial variables, the corresponding correlation coefficients result from the statistical properties of the gridded rainfall depth quantified using 1,000 simulations of rainfall depths at each grid (see Figure 18). From Figure 18, the correlation coefficient of the mean value and the upper bound of 95% confidence interval exceeds 0.8; this means that the averages of the simulated gridded rainfall depths are strongly related to the observations in space. Also, in spite of the larger relative errors regarding the standard deviation and lower bound of 95% confidence interval calculated from the simulations, about −1.47 and 5.04, their spatial correlation coefficients approximate 0.5.

To sum up the above results, the average rainfall depths at various grids have good agreement with those derived from observations. However, their standard deviations are mostly less than those from the observations which leads to a narrow 95% confidence interval consisting of the higher lower bound and smaller upper bounds, respectively, probably in order to persist in the correlation structures regarding the rainfall depths in space. In detail, as for the upper bound of the 95% confidence interval, the maximum (about 770 mm) can be found at the 336th grid and it slightly exceeds the observation (approximately 680 mm), implying that the proposed SM_STGR model is advantageous in the generation of heavy rainstorms probably applied in the investigations related to the flooding and inundation simulation. Consequently, the proposed SM_GSTR model can reproduce the simulations of event-based gridded rainfall depths with a high spatial correlation in comparison to those from the observations.

#### Gridded storm patterns

Similar to the gridded rainfall depths, the performance evaluation of simulating gridded storm patterns is accomplished by comparing the statistical properties and corresponding dependence in time and space between observation and simulations. Since the simulation of the gridded storm pattern is combined with the generated areal average of cumulative dimensionless rainfalls and associated deviation at dimensionless times, the performance assessment of simulating the gridded rainstorm can be achieved in comparison to the statistical properties of the areal average and associated gridded deviation, respectively. The relevant results are addressed and assessed below.

##### Areal average

Figure 19 shows the comparison of the statistical properties of the areal average of the cumulative dimensionless rainfalls from the simulations with those from the observation; the statistical properties of the simulations excellently match those from the cumulative dimensionless rainfall with the low relative errors being −0.095 (mean), −0.263 (standard deviation), 0.76 (upper bound), and −0.078 (lower bound). These results reveal that the proposed SM_GSTR model can reproduce the areal average of storm pattern with high reliability by permitting well the temporal statistical properties of the storm pattern.

Since the areal average of dimensionless rainfalls can be treated as the temporal variate, the correlation coefficients of their statistical properties can be used to assess the dependence between the observed and simulated areal averaged dimensionless rainfalls. Figure 20 presents the comparison in the correlation coefficients of the statistical properties, indicating that the above correlation coefficients mostly reach 0.95. This result reveals that the varying trend of the simulated areal average of dimensionless rainfalls in time exhibits a good match with the observations. This matching indicates that the proposed SM_GSTR model can generate a large number of areal averaged storm patterns of which the temporal dependence is consistent with the observations.

##### Gridded deviation

The assessment of generating the gridded deviations of the cumulative dimensionless rainfalls at the dimensionless times of are made as shown in Figure 21. It can be known that the mean and standard deviation of the simulated grid-based bias at the four dimensionless times have an acceptable agreement with those from observations in which the relative errors are on average −0.011 and −0.045, respectively. Moreover, the fitness of the upper bound of 95% confidence interval to the results from observations is better than the lower bound due to the lower bound of 95% confidence interval simulated gridded bias being mostly overestimated. However, it can be seen that the spatial change in the bounds of 95% confidence interval has good agreement with the observation as a result of the high correlation coefficient being 0.6.

In addition to the areal average of dimensionless rainfalls, this study quantifies the correlations among the gridded deviations of the cumulative dimensionless rainfall at the four dimensionless times as shown in Figure 22. Figure 22 shows that the mean value of the simulated gridded deviations at the four dimensionless times is evidently and strongly related to that from observation with the correlation coefficient being 0.9. As a result, even for the high statistical moments (the standard deviation and 95% confidence interval) (see Figure 22), there is good fitness to the results from the observed values of the gridded deviation, because the correlation coefficients between the simulated and observed deviation exceed 0.6.

##### Combination of areal average and gridded deviation

Furthermore, according to a graphical comparison of storm patterns at the specific grids for various simulated grid-based rainstorm events (see Figure 23), there are similar temporal changes in storm patterns at the six grids randomly selected in the Nankan River watershed. In detail, the six grids concerned exhibit similar types of storm patterns, with different gridded deviations. Especially for the neighboring grids (e.g., Grid185 and Grid188), the temporal varying trend of storm pattern at Grid185 is comparable to that at Grid188, and their corresponding difference is, on average, merely about 0.0016.

To sum up the above results, the proposed SM_GSTR method can take into account the statistical properties of gridded rainstorms, especially for the spatiotemporal correlation, to generate the time series of rainfall at all grids within a study area. Also, the resulting simulations of the gridded rainfall characteristics by the proposed SM_GSTR method are markedly consistent with the results from the observations, except for the standard deviation of the gridded rainfall depth and the upper bound of the 95% confidence interval. As a result, the proposed SM_GSTR method can well capture the behavior in time and space of gridded rainstorms.

### Spatial comparison of simulations of gridded rainstorms

Eventually, the generated grid-based rainfall characteristics can be combined as the simulated hyetographs at all grids. Figure 24 illustrates the 1st, 100th, and 1,000th simulated rainstorms at 336 grids in the Nankan River watershed, comprising the generated gridded rainfall characteristics. It can be found that the above four gridded rainstorm events are significantly associated with different changes in space and time. Also, in the case of spatial distribution of hourly gridded rainfall for the 13-hour simulated event as shown in Figure 25, the rainfall reaches the maximum around the 11th hour and the grids with larger rainfall are mostly located near the lower right corner of the study area. It summarizes that the spatial change in the gridded rainstorm simulated by the proposed SM_GBTR model varies with time and resembles the spatial distribution of real rainstorms.

## CONCLUSION

This study aims to present a stochastic method (SM_GSTR) for generating a large number of gridded short-time (hourly) rainstorm events within a watershed/urban area by means of the correlated MMCS approach (Wu *et al.* 2006). In the proposed SM_GSTR model, the rainstorms at all grids of interest consist of the four spatially and temporally correlated rainfall characteristics, the rainfall duration, gridded rainfall depths, and gridded storm pattern composed of the areal average storm patterns and the associated gridded deviation of storm pattern.

A basin located in north Taiwan, the Nankan River watershed, is selected as the study area and the associated grid-based precipitation data provided by Center Weather Bureau in Taiwan of 20 typhoon events are used to demonstrate its feasibility and applicability, in which the statistical properties of gridded rainfall characteristics are compared to those from the observations. The results from the model feasibility indicate that for the comparison of the statistical properties of gridded rainfall characteristics, the mean values of gridded rainfall depths and the mean and standard deviation of the areal average storm pattern have excellent agreement with those from their observation. Also, the correlation coefficients of the mean value and upper bound of 95% confidence interval of grid-based rainfall characteristics (about 0.8) are quite consistent with those from observations, whereas the correlation coefficient of the standard deviation and lower bound of 95% confidence interval (approximately 0.65) merely differ from observations attributed to an acceptable bias. It is expected that the proposed SM_GSTR can reproduce a large number of the gridded rainfall characteristics, of which the corresponding statistical features approximate those from observations with an acceptable and reasonable bias.

The proposed SM_GSTR model is established by modifying Wu's method under consideration of the temporal and spatial dependence structure of rainstorms with various temporal and spatial scales. However, the effect of short-range scale regarding the fractal behavior of rainfall should impact the variation in the precipitation and induced runoff (e.g., Kantelhardt *et al.* 2006). Therefore, the effect of multifractal temporal scale on the simulation of storm pattern would be evaluated in future work. Also, due to the overestimated lower bound of the 95% confidence interval regarding the gridded rainfall depths, the proposed SM_GSTR model might have difficulty generating the gridded rainfall depths of light rainstorms. Therefore, further work should be done by improving the SM_GSTR model under consideration of persistence of the standard deviations regarding the gridded rainfall depths.

In addition to improvement regarding the model development, it is well known that climate change impacts the evaluation in the hydrological analysis so strongly as to require more reliable estimates of future precipitation (Olsson *et al.* 2013). Nevertheless, the proposed SM_GSTR method reproduces a large number of high-resolution rain fields comprising the gridded rainstorm events based on the statistical properties of observed QPE data without considering the effect of the climatic conditions (e.g., temperature, moisture, and CO_{2}) in response to climate change (Doroszkiewicz *et al.* 2019). Therefore, by integrating with a physically based climatic simulation model, such as the GCMs, the proposed SM_GSTR model would be capable of producing the near-real simulations of the gridded rainstorm events under the consideration of climate change (e.g., Qiang *et al.* 2020). Furthermore, the proposed SM_GSTR method would be combined with the relevant hydrologic/hydraulic numerical simulation models to reproduce massive rainfall-induced simulation scenarios with high resolution in time and space (e.g., flood and inundation); the resulting gridded rainfall-induced simulations could be applied in training an AI model for flood forecasts.

## DATA AVAILABILITY STATEMENT

All relevant data are included in the paper or its Supplementary Information.

## REFERENCES

_{i}=

*γ*+

*δ*ln[y

_{i}

^{α}/(1-y

_{i}

^{α})] by methods of translation