10.1 Introduction

Rainfall frequency analysis is of fundamental importance for hydrologic and hydraulic engineering design. In what follows, we will first introduce some examples with regard to rainfall analysis. Rainfall intensity-duration-frequency (IDF) or rainfall depth-duration frequency (DDF) curves published by National Oceanographic Atmospheric Administration (NOAA) are classic examples of rainfall frequency analysis. The IDF (or DDF) curves are derived first by separating rainfall events based on their durations (e.g., 15 minutes, 30 minutes, one hour, etc.) and then by fitting a univariate probability distribution to the rainfall depth or intensity data of a certain duration. The fitted univariate distribution is applied to produce a family of rainfall depth-frequency curves. In this manner, the two-dimensional depth-duration analysis is reduced to a one-dimensional analysis, involving only intensity (or depth) corresponding to a fixed duration. As described by the NOAA documents (e.g., TP-40), the IDF (or DDF) curves may be estimated from either annual maximum series or partial duration series. The IDF (or DDF) curves are widely applied in hydrological and hydraulic engineering design.

The rational method relates rainfall intensity (I) of a given duration (normally equal to the time of concentration) of a certain return period to peak runoff (discharge) (Q), where the peak runoff is assumed as a linear function of rainfall (Q = CIA)Q=CIA), where A is the area of the drainage basin. In this method, rainfall of a certain return period results in the runoff peak of exactly the same return period. To date, the rational method is commonly applied in urban hydrology (e.g., urban rainfall and runoff analysis) and urban hydraulic engineering design (e.g., detention/retention basin design, storm sewer design, and highway drainage design).

The SCS method, developed by Soil Conservation Service (now, the Natural Resources Conservation Service), may be applied to larger areas compared to the rational method (usually less than 60 acres [about 25 hectares]) for estimating runoff of a given rainfall amount. This method estimates the amount of surface runoff (or excess rainfall) through what is called the Curve Number (CN), which is related to land use and land cover, antecedent soil moisture, hydrologic condition, and soil moisture retention capacity.

The probable maximum precipitation (PMP) method, which does not rely on the IDF (DDF) curve, estimates the maximum amount of precipitation that may probably occur. The PMP analysis is required for the design of dams, dam breach analysis, spillway analysis, design of nuclear power plants, etc. These examples may be considered to illustrate applications of univariate rainfall analysis in hydrologic and hydraulic engineering design.

In the past three decades, bivariate (and multivariate) rainfall frequency analysis has attracted significant attention, because rainfall variables may be correlated and may significantly affect surface runoff (Cordova and Rodriguez-Iturbe, 1985). In the early days, the bivariate exponential distribution was applied to model the correlation structure of extreme rainfall variables (e.g., Hashino, 1985; Singh and Singh, 1991; Bacchi et al., 1994). Later, other bivariate rainfall models were investigated to model the relation between rainfall intensity and rainfall duration, for example, improved derived flood frequency distribution (DFFD) model by Kurothe et al. (1997) and Goel et al. (2000); Yue (2000a, 2000b, 2000c) investigated the applicability of bivariate normal, Gumbel logistic, and Gumbel mixed distributions. Besides the application to river discharge (Favre et al., 2004), the copula theory has been applied to bivariate and multivariate rainfall analysis (Grimaldi et al., 2005; Zhang and Singh, 2007a, 2007b, 2007c; Kao and Govindaraju, 2007, 2008; Cong and Brady, 2012; Zhang et al., 2012; Hao and Singh, 2013; Zhang et al., 2013; Abdul Rauf and Zeephongsekul, 2014; Cantet and Arnaud, 2014; Khedun et al., 2014; Moazami et al., 2014; Vernieuwe et al., 2015; among others).

With the advantages of the copula theory discussed in the preceding chapters, we will illustrate the application of copula theory to bivariate (or multivariate) rainfall frequency analysis. It is assumed that rainfall variables are continuous variates. However, rainfall variables may actually be discrete in nature.

10.2 Rainfall Depth-Duration Frequency (DDF) Analysis

Many studies have employed copulas for bivariate (multivariate) rainfall analysis based on annual maximum series (AMS). In this section, we will use Partial Duration Series (PDS) to illustrate the copula application to derive the DDF curves. The rainfall data with a 15-minute interval were collected for the rain gauge station: coop-166394 near Morgan City, Louisiana. The recorded data cover a period from May 8, 1971, to January 1, 2014. The rainfall data are available upon request from National Climate Data Center (NCDC).

The general procedure for DDF analysis includes the following steps:

1. 
   

   1. Separate the rainfall records collected into independent rainfall events. Extract the rainfall depth and rainfall duration from these independent rainfall events obtained.

   
   
2. 
   

   2. Evaluate the marginal rainfall depth and rainfall duration variables and corresponding marginal distributions.

   
   
3. 
   

   3. Evaluate the rank-based correlation of rainfall depth and rainfall duration. Choose the possible copula candidates.

   
   
4. 
   

   4. Perform the rainfall depth and rainfall duration analysis with the use of the possible copula candidates. Select the best-fitted copula functions.

   
   
5. 
   

   5. Estimate the rainfall depth of given rainfall duration for a given return period.

• 
  

  In what follows, we will discuss how to perform the DDF analysis in detail.

10.2.1 Rainfall Data Processing

Before analyzing bivariate rainfall variables (i.e., rainfall depth and duration), we need to separate the rainfall data into individual rainfall events first. As commonly done, a six-hour duration of no rain is considered as the criterion to separate any two events. From a total of 12,089 available rainfall records for rain gage coop-166394, a total of 2,816 events were identified for the 43-year duration. Table 10.1 illustrates the rainfall event separation year 1971 as an example. From Table 10.1, it can be seen that there are nine independent rainfall events identified from May 22, 1971, to June 29, 1971. Table 10.2 lists the nine rainfall events separated. As an example, we consider the No. 5 event, which started on June 20, 1971, at 14:15 and ended on 15:15 on the same day. Summing up the incremental rainfall depths within this time window, we have the following:

depth=2.54+2.54+10.16+7.62=22.86mm.maximum rainfall intensity=10.16/0.25=40.64mm/h.duration=1h.

Table 10.1. Illustration of rainfall event separation.

Event no.	Date	Rainfall amount (mm)a	Interarrival time (h)b
	19710522 15:30	7.62	—
	19710601 00:15	0	—
1	19710605 16:15	2.54	336.75
2	19710616 13:45	2.54	261.50
3	19710618 16:00	5.08	50.25
	19710618 17:30	2.54	1.50
4	19710619 13:30	2.54	20.00
	19710619 13:45	2.54	0.25
	19710619 14:30	2.54	0.75
5	19710620 14:30	2.54	24.00
	19710620 14:45	2.54	0.25
	19710620 15:00	10.16	0.25
	19710620 15:15	7.62	0.25
6	19710621 10:00	2.54	18.75
	19710621 10:15	7.62	0.25
7	19710622 20:00	12.7	33.75
	19710622 20:15	5.08	0.25
	19710622 20:30	5.08	0.25
	19710622 20:45	5.08	0.25
	19710622 21:00	2.54	0.25
	19710622 21:15	2.54	0.25
	19710622 22:00	2.54	0.75
	19710622 22:30	2.54	0.50
8	19710624 15:30	2.54	41.00
9	19710629 12:00	5.08	116.50

Notes:

^a Incremental rainfall depth with 15-minute interval until the time stated.

^b Difference between current day and time with the previous day and time.

Table 10.2. Nine rainfall events separated based on Table 10.1.

No	Depth (mm)	Duration (h)	Max. Intensity (mm/h)a	Start	End
1	2.54	0.25	10.16	6/5/71 16:00	6/5/71 16:15
2	2.54	0.25	10.16	6/16/71 13:30	6/16/71 13:45
3	7.62	1.75	20.32	6/18/71 15:45	6/18/71 17:30
4	7.62	1.25	10.16	6/19/71 13:15	6/19/71 14:30
5	22.86	1.00	40.64	6/20/71 14:15	6/20/71 15:15
6	10.16	0.50	30.48	6/21/71 9:45	6/21/71 10:15
7	38.1	2.75	50.8	6/22/71 19:45	6/22/71 22:30
8	2.54	0.25	10.16	6/24/71 15:14	6/24/71 15:30
9	5.08	0.25	20.32	6/29/71 11:45	6/29/71 12:00

Note: ^a Maximum average intensity of 15-minute interval.

Similarly, all rainfall events may be separated based on the six-hour duration of no rain criterion. Setting the threshold for the identified rainfall events as follows:

Threshold = median (rainfall depth) + std (rainfall depth)

Threshold=medianrainfall depth+stdrainfall depth(10.1)

From the record, we have median = 7.62 mm and standard deviation = 23.87 mm, which yield the threshold = 31.49 mm. Applying this threshold, we reduced the number of rainfall events to 378 that is roughly about nine events per year. With the partial duration rainfall series thus identified, we can then start to investigate bivariate rainfall characteristics through (i) the investigation of the marginal distribution and (ii) the investigation of dependence.

10.2.2 Investigation of Marginal Distributions: Depth and Duration

Before we investigate marginal distributions, we will first look at a scatter plot of rainfall variables (Figure 10.1(a)). Zooming in on the lower-right corner (Figure 10.1(b)), we see there are ties in both rainfall depth and rainfall duration variables. The kernel density (nonparametric probability density estimation (Wand and Jones, 1995) is applied to approximate the nonparametric probability density and distribution function for univariate rainfall variables. The kernel density function is given as follows:

f̂x=1nh∑i=1nKx−xih(10.2)

In Equation (10.2), K(.)K. is the kernel function. Here we use the commonly applied K(x) = ϕ(x)Kx=ϕx, i.e., the normal kernel (the normal density function); h is the smoothing parameter, which is also called bandwidth (h = 6.086mm, 1.797hrh=6.086mm,1.797hr for rainfall depth and rainfall duration respectively); and n is the sample size.

Figure 10.1 Scatter plot for rainfall depth and rainfall duration: (a) original; (b) zoomed in at lower-right corner.

To compute the probability density and marginal probability using the kernel density, the MATLAB function is applied as follows:

pdf = ksdensity(x, x₁,^‘support^‘,^‘positive^‘)

pdf=ksdensityxx1,′support′,′positive′(10.2a)

cdf = ksdensity(x, x₁,^‘function^‘,^‘cdf^‘,^‘support^‘, ^‘positive^‘)

cdf=ksdensityxx1,′function′,′cdf′,′support′′positive′(10.2b)

In Equations (10.2a) and (10.2b), x and x₁ represent the random variable and the data points where the nonparametric pdf and cdf need to be evaluated. Figure 10.2 plots the density function as well as the cumulative probabilities for both rainfall variables. The CDF estimated from the kernel density is applied for bivariate analysis using copulas.

Figure 10.2 Frequency and cumulative probability plots with kernel density function for rainfall depth and rainfall duration series.

10.2.3 Bivariate Rainfall Frequency Analysis

The scatter plot in Figure 10.1 indicates positive dependence between rainfall depth and rainfall duration. The rank-based Kendall correlation coefficient is computed as τ_n≈0.32τn≈0.32. Among the copula candidates (i.e., Gumbel–Hougaard, Clayton, Frank, Gaussian, and Student’s t copulas), the Frank copula is found to better model the bivariate rainfall characteristics. Figure 10.3 compares the empirical CDF estimated from the kernel density with the bivariate random variables simulated from the fitted Frank copula with its parameter value of 3.529. Comparison shows that (i) the simulated random variates cover the overall dependence fairly well; and (ii) the tie existing in both rainfall depth and rainfall duration variables may impact the concordance of the bivariate rainfall variables.

Figure 10.3 Comparison of bivariate empirical distribution using kernel density with the random variables simulated from the fitted Frank copula.

However, with the continuous assumption, we will proceed to estimate the rainfall depth for a given duration of a given return period. The exceedance probability (P_exPex) corresponding to a given return period (T) for the partial duration series may be written as follows:

Pex=1μT(10.3)

In Equation (10.3), μ≈9μ≈9, the average number of events per year.

Equating Equation (10.3) to the exceedance probability of rainfall depth of a given rainfall duration, we have the following:

Praindep>xRaindur=d=1μT(10.4)

Equation (10.4) is equivalent to the following:

CFrankFdep≤FdepxFdur=Fdurd=1−1μT(10.5)

In Equation (10.5), C^Frank(F_dep ≤ F_dep(x)| F_dur = F_dur(d)) = P(dep ≤ x| dur = d)CFrankFdep≤FdepxFdur=Fdurd=Pdep≤xdur=d. The conditional copula in Equation (10.5) is listed as #5 in Table 4.2. Applying the kernel density to the given durations of 1, 2, 3, 6, 12, and 24 fours, we have F_dur(d)Fdurd computed as follows:

F_dur(d) = [0.0818, 0.1385, 0.2079, 0.4362, 0.7480, 0.9551]

Fdurd=0.08180.13850.20790.43620.74800.9551.

For the return period of 1, 2, 5, 10, 25, 50, and 100 years, we have the exceedance probability computed using Equation (10.3) directly as follows:

P_ex = [0.8862, 0.9431, 0.9772, 0.9886, 0.9954, 0.9997, 0.9989]

Pex=0.88620.94310.97720.98860.99540.99970.9989.

Substituting F_dur(d), P_exFdurd,Pex into Equation (10.5), we can compute F_dep(x)Fdepx numerically using the bisection method. Finally, we can estimate the corresponding rainfall depth using the inverse of the kernel density (fitted to the observed rainfall depth) with the computed F_dep(x)Fdepx. Table 10.3 lists the estimated F_dep(x)Fdepx and the corresponding estimated rainfall depth. Figure 10.4 compares the rainfall depth estimated from copula-based analysis with the published DDF of partial duration for Morgan City, Louisiana (http://hdsc.nws.noaa.gov/hdsc/pfds/pfds_map_cont.html?bkmrk=la). Comparison shows that (i) for the storms with shorter duration and return periods less than 10 years, the copula estimates are either closely following the NOAA estimates or well within NOAA 90% bounds; (ii) for short durations (i.e., D = 1 and 2 hours) and higher return periods (T ≥ 25 yr)T≥25yr), the copula estimates are higher than the NOAA upper 90% bounds; and (iii) as the storm duration increases, the copula estimates for higher return periods get closer to either NOAA upper 90% bounds or actually closely follow the NOAA estimates.

Table 10.3. Estimated probability distribution of rainfall depth and estimated rainfall depth of given duration with given return period.

	1-yr	2-yr	5-yr	10-yr	25-yr	50-yr	100-yr
				Fdep(x)Fdepx
1-hr	0.5988	0.7364	0.8656	0.9252	0.9677	0.9834	0.9916
2-hr	0.6512	0.7793	0.8921	0.9413	0.9751	0.9873	0.9936
3-hr	0.7034	0.8195	0.9153	0.9548	0.9811	0.9904	0.9952
6-hr	0.8296	0.9064	0.9599	0.9794	0.9916	0.9958	0.9979
12-hr	0.9250	0.9622	0.9848	0.9924	0.9969	0.9985	0.9992
24-hr	0.9594	0.9802	0.9922	0.9961	0.9984	0.9992	0.9996
	Rainfall depth (mm)
1-hr	55.30	68.64	92.78	116.32	166.01	206.38	246.33
2-hr	59.54	74.73	101.12	128.39	182.26	221.99	261.82
3-hr	64.72	81.89	110.86	143.97	198.76	238.62	277.24
6-hr	84.01	106.74	151.78	193.73	246.69	284.32	314.98
12-hr	116.25	155.60	211.31	251.90	299.25	327.08	350.48
24-hr	150.98	195.86	250.48	287.99	326.54	350.01	370.53

Figure 10.4 Comparison of copula estimates with the NOAA estimations with a 90% confidence bound.

The differences between the NOAA-DDF and the copula-based DDF curves may be due to the following:

1. 
   

   i. The NOAA-DDF analysis only extracts rainfall events for certain durations. These extracted events are then treated as univariate random variables and are fitted by univariate probability distributions.

   
   
2. 
   

   ii. In the copula-based DDF analysis, on the other hand, rainfall events extracted may yield different rainfall durations. The bivariate rainfall depth-duration model is then constructed, and the rainfall depth of a given duration is estimated from the conditional probability function of f(depth < depth^∗| duration = duration^∗)fdepth<depth∗duration=duration∗. In this analysis, the duration can take on any value.

   
   
3. 
   

   iii. The ties that may exist in the NOAA-DDF extracted events may not have the same degree of impact as that of copula-based DDF events. As discussed earlier, there may be many ties in the rainfall depth and duration of the extracted rainfall events (partial duration or annual maximum series), and these tied values may distort the concordance of the bivariate rainfall analysis. Additionally, the rainfall variables (especially rainfall duration) may be discrete in nature.

Even with the differences between the NOAA and copula-based DDF curves constructed for the partial duration time series, the copula-based method may be considered as a rational alternative for rainfall DDF (or IDF) construction with simpler and faster rainfall separation (events regardless of the length of rainfall duration) compared to that of NOAA analysis (rainfall duration–based directly).

10.3 Spatial Analysis of Annual Precipitation

With the assumption of annual precipitation amount as a random variable, the general procedure for spatial analysis of annual precipitation includes the following steps:

1. 
   

   1. Select the region of interest, identify the rain gauges, and collect the annual precipitation records.

   
   
2. 
   

   2. Evaluate the pairwise rank-based correlation coefficient of annual precipitation.

   
   
3. 
   

   3. Identify the possible vine structure based on the rank-based correlation coefficients computed, and select possible copula candidates for T1 first, and then proceed with the analysis for the rest of the tree structure as discussed in Chapter 5.

   
   
4. 
   

   4. Identify the proper tree structure for the asymmetric Archimedean copula and then proceed with the analysis as discussed in Chapter 5 for the asymmetric Archimedean copula.

   
   
5. 
   

   5. Construct the meta-elliptical copula for the multivariate precipitation variables.

   
   
6. 
   

   6. Compare the performance of different copula construction approaches.

To illustrate the spatial analysis of annual precipitation (rainfall), we will use four NOAA rainfall stations located in the Cuyahoga River Watershed, Ohio (see Table 10.4). The copula model is constructed from the annual rainfall data collected from 1953 to 2012 from NCDC. In this case study, we will apply D-vine, meta-elliptical copulas (i.e., meta-Gaussian and meta-Student T) and asymmetric Archimedean copulas. The reason that a D-vine copula is chosen from the pair copula construction is that there is no obvious center variable governing the dependence structure among all four rainfall stations (see the rank-based Kendall correlation coefficient listed in Table 10.5).

10.3.1 Application of D-Vine Copula to Four-Dimensional Rainfall Variables

Copula Identification for T₁

According to Kendall’s tau correlation coefficient matrix, the proper structure for T₁ is as follows: R330058 − R333780 − R336949 − R331458R330058−R333780−R336949−R331458 (i.e., the bivariate pairs for T₁ are [R330058, R333780]; [R333780, R336949]; [R336949, R331458]). Using the empirical marginals (Weibull plotting position formula), let U₁, U₂, U₃, U₄U1,U2,U3,U4 represent the empirical marginals as follows:

U1=F̂nR330058;U2=F̂nR333780;U3=F̂nR336949;andU4=F̂nR331458.

The D-vine structure for this example is the same as in Figure 10.5. In this case study, we choose Archimedean copulas for dealing with the positive dependence (Gumbel–Hougaard, Clayton, Frank, Joe, and BB1 copulas) as the candidates. Chapter 4 listed the one-parameter Archimedean copulas candidates. Hence we only give the formula for BB1 copula, which is a two-parameter Archimedean copula with the limiting conditions of either the Clayton or Gumbel–Hougaard copula. The BB1 copula (Joe, 1997) can be formulated as follows:

Cuvθ1θ2=1+u−θ1−1θ2+v−θ1−1θ21θ2−1θ1;θ1>0;θ2≥1(10.6)

The BB1 copula converges to (i) the Gumbel–Hougaard copula if θ₁→0θ1→0; and (ii) the Clayton copula if θ₂ = 1θ2=1.

Figure 10.5 D-vine structure for four-dimensional rainfall variables: (1) R330058, (2) R336949, (3) R333780, and (4) R331458.

In addition, the BB1 copula has both upper- and lower-tail dependence coefficients, as follows:

λL=2−1θ1θ2,λU=2−21θ2

The parameters of T₁ are estimated with the pseudo-MLE through the empirical marginals for all the copula candidates (Table 10.5). Table 10.6 also lists the log-likelihood, AIC, and BIC values with the best-fitted copula highlighted. From Table 10.6, we see that the two-parameter BB1 copula is the best-fitted copula for stations (R330058, R333780, R333780, and R336949), and the Gumbel–Hougaard copula is the best-fitted copula for stations R336949 and R331458.

Table 10.6 Estimation results for copula candidates.

Variables	Copulas	θθ	LL	AIC	BIC
	Gumbel-Hougaard (GH)	2.7782	35.0603	−68.1206	−66.0601
	Clayton (C)	2.9226	33.9990	−65.9979	−63.9375
U1 v. s. U2U1v.s.U2	Frank (F)	8.8627	31.7661	−61.5322	−59.4718
	Joe (J)	3.3021	28.8828	−55.7655	−53.7051
	BB1	[1.0203, 1.9788]	39.0411	−74.0821	−69.9612
	Gumbel-Hougaard (GH)	2.3336	26.0645	−50.1290	−48.0685
	Clayton (C)	1.8549	20.2556	−38.5113	−36.4508
U2 v. s. U3U2v.s.U3.	Frank (F)	6.8861	22.8811	−43.7622	−41.7017
	Joe (J)	2.8627	23.0147	−44.0294	−41.9689
	BB1	[0.3841,2.0196]	26.7981	−49.5963	−45.4754
	Gumbel-Hougaard (GH)	1.7527	22.9526	−43.9052	−41.8447
	Clayton (C)	1.2274	12.7401	−23.4801	−21.4197
U3 v. s. U4U3v.s.U4	Frank (F)	4.9432	14.3614	−26.7228	−24.6623
	Joe (J)	1.9670	10.9197	−19.8394	−17.7790
	BB1	[0.4963,1.4682]	15.1711	−26.3423	−22.2214

Only gold members can continue reading. Log In or Register to continue

Share this:

• Click to share on Twitter (Opens in new window)
• Click to share on Facebook (Opens in new window)

Related

Related posts:

1. 3 – Copulas and Their Properties
2. 7 – Non-Archimedean Copulas
3. 13 – Drought Analysis
4. 12 – Water Quality Analysis