Quantifying uncertainty in urban flooding analysis considering hydro-climatic projection and urban development effects

How will the combined impacts of land use change, climate change, and hydrologic modeling influence changes in urban flood frequency and what is the main uncertainty source of the results? Will such changes differ by catchment with different degrees of current and future urban development? We attempt to answer these questions in two catchments with different degrees of urbanization, the Fanno catchment with 84% urban land use and the Johnson catchment with 36% urban land use, both located in the Pacific Northwest of the US. Five uncertainty sources – general circulation model (GCM) structures, future greenhouse gas (GHG) emission scenarios, land use change scenarios, natural variability, and hydrologic model parameters – are considered to compare the relative source of uncertainty in flood frequency projections. Two land use change scenarios, conservation and development, representing possible future land use changes are used for analysis. Results show the highest increase in flood frequency under the combination of medium high GHG emission (A1B) and development scenarios, and the lowest increase under the combination of low GHG emission (B1) and conservation scenarios. Although the combined impact is more significant to flood frequency change than individual scenarios, it does not linearly increase flood frequency. Changes in flood frequency are more sensitive to climate change than land use change in the two catchments for 2050s (2040–2069). Shorter term flood frequency change, 2 and 5 year floods, is highly affected by GCM structure, while longer term flood frequency change above 25 year floods is dominated by natural variability. Projected flood frequency changes more significantly in Johnson creek than Fanno creek. This result indicates that, under expected Correspondence to: H. Chang (changh@pdx.edu) climate change conditions, adaptive urban planning based on the conservation scenario could be more effective in less developed Johnson catchment than in the already developed Fanno catchment.


618
I.-W. Jung et al.: Quantifying uncertainty in urban flooding analysis considering hydro-climatic projection large uncertainty, identifying potential changes in flood risk according to climate and land use changes is an important area of concern to water resource managers and land use planners (Hine and Hall, 2010).
For mitigation of and protection from potential flood risk in urban areas, we need to improve our understanding of the possible impacts of the ubiquitous uncertainty of urban flood projection. This uncertainty stems from several sources; internal variability of the climate system, future GHG and aerosol emissions, the translation of these emissions into climate change by GCMs, spatial and temporal downscaling, and hydrological modeling . Uncertainty will not be radically removed or reduced until the development of the new technology of climate and hydrologic modeling based on additional observation of hydrometeorological variables, such as soil moisture, snow, actual evapotranspiration, and groundwater. This uncertainty complicates the accurate interpretation of climate impact assessment. Therefore, many researchers have attempted to quantify the irreducible uncertainty in hydrologic streamflow projections (New et al., 2007;Wilby, 2005;Kingston and Taylor, 2010), low flow (Wilby and Harris, 2006), flooding (Booij, 2005;Kay et al., 2009;Raff et al., 2009;Moradkhani et al., 2010), and drought (Ghosh and Mujumdar, 2007;Mishra and Singh, 2009). Despite substantial effort of previous studies, however, large uncertainty in climate impact studies still remain .
Floods in urban areas are controlled by the integrated condition of geophysical characteristics, urban infrastructure, drainage system, and hydro-climatologic regime (Epting et al., 2009). Thus, different levels of urban development could lead to different hydrologic responses among catchments, though they are under identical climate change ). Kay et al. (2009) investigated the uncertainty in climate change impact on flood frequency for two catchments in England, showing that uncertainty can vary significantly between catchments that have different rainfall regimes and topographic characteristics. Prudhomme and Davies (2009) reported similar findings for four catchments in Britain. Additionally, the combined effects of climate change and anthropogenic land use change significantly aggravate the accuracy of hydrologic prediction associated with overall urban environmental management (Brath et al., 2006;Choi, 2008;Franczyk and Chang, 2009;Praskievicz and Chang 2009a;Tu, 2009). However, relatively few studies have examined the combined effects of climate change and urban development on the uncertainty of urban flood projections in catchments with different degrees of urban development. This study attempts to fill this gap using two future land use change scenarios projected under two GHG emission scenarios.
The three research questions are: (1) What are the main sources of uncertainties affecting the changes in urban flood frequency? (2) How will the combined impacts of land use change and climate change influence changes in flood frequency? and (3) How is flood frequency projected to change in two urban catchments with different degrees of urban development for the 2050s (2040-2069) with respect to the reference period 1960-1989? This paper can contribute to a better understanding of the combined impact of climate and land use changes on urban flood frequency, and is expected to help decision makers with practical urban planning and management to mitigate potential flood damage in urban areas in a changing climate.

A process of flood frequency uncertainty analysis
We investigate changes in flood frequency and the uncertainties associated with the combined effects of climate change and land use change in two catchments -Fanno Creek (80.5 km 2 ) and Upper Johnson (hereafter Johnson) Creek (68.3 km 2 ) in the Portland metropolitan area of Oregon, USA. The Fanno catchment is highly developed with 84% urban land use, and the Johnson catchment is moderately developed with 36% urban land use in 2001 (see Fig. 1).
To quantify uncertainty in flood frequency change, this study considers five uncertainty sources; GCM structures, future GHG emission scenarios, future land use scenarios, hydrologic model parameters, and natural variability of the climate system. The GCM simulations are downscaled using the delta method to correct the bias between simulated and observed precipitation and temperature, which is attributed to scale mismatch between GCMs and catchment hydrologic models, as well as the lack of sub-grid scale climate dynamics such as orographically convective precipitation (Im et al., 2010). Precipitation Runoff Modeling System (PRMS), a physically-based, deterministic, and semi-distributed model, is employed to simulate daily runoff changes and resulting changes in flood frequency under different climate and land use conditions. PRMS has been applied successfully in several regions with varying climate and land use (Bae et al., 2008a;Clark et al., 2008;Hay et al., 2006;Qi et al., 2009;Viney et al., 2009). To better understand the wide array of individual and combined factors that can affect the hydrologic response in a watershed system, Risley et al. (2010) employed PRMS, driven by GCM outputs, in 14 watersheds across the US, and conducted a comparative statistical analysis on the outputs. In the Willamette River basin, Oregon, PRMS has been applied in a water quality study (Laenen and Risley, 1997) and in a climate change impact study . To consider PRMS model parameter uncertainty, we extract acceptable parameter sets based on the Nash-Sutcliffe efficiency (NSE) criterion that estimates the degree of closeness between observed and simulated streamflow. Latin Hypercube Sampling is employed to efficiently sample the PRMS parameter sets within plausible ranges. A similar approach was undertaken by Wilby and Harris (2006).
It is also important to find whether the changes in flood frequency for the future period are larger than the natural (or model internal) climate variability (Hagemann and Jacob, 2007). It is especially likely that precipitation change derived from different initial conditions of GCMs could lead to different interpretation of the results due to large natural internal variability. To estimate natural climate variability, we employ the moving block Bootstrap resampling method (Ebtehaj et al., 2010), which produces a large number of new climate series through random selection of observed climate data. This method allows us to explore the range of different flood frequencies that could be obtained by our finite sampling of the internal climate variability (Kay et al., 2009). The US Geological Survey's PeakFQ program (Flynn et al., 2006) is applied to estimate flood frequency with different recurrence intervals such as 2, 5, 10, 25, 50, and 100 years. To represent realistic future land use changes, we use two land use change scenarios: the conservation and the development scenarios, developed by the PNWERC (2002), these are both compared with 2001 land use. Further details of the data and methods used in this study are given in the following sections.

Study area and data
Fanno creek and the Johnson creek are important resources in the Portland metropolitan area, located in the valley of the Willamette River basin in Oregon (see Fig. 1). As a source of recreation and wildlife (Laenen and Risley, 1997), they contribute to the regional socio-economic and environmental systems. Two catchments are located in a modified marine temperate climate region in which summers are warm and dry but winters are cold and wet. More than 80% of the annual precipitation occurs from October through May and less than 10% precipitation falls during July and August (Praskievicz and Chang, 2009b). This seasonality of precipitation causes periodic flooding and companying travel disruptions in winter .
In our study areas, most precipitation is in the form of rainfall. Unusual snow melts quickly during subsequent rain storms (Lee and Snyder, 2009). Therefore, the surface hydrology of these regions is highly dominated by frequent rainfall. Although Fanno and Johnson are close to each other and have identical climate conditions, they show different hydrologic regimes. Fanno shows a higher runoff ratio, defined as the ratio of total monthly runoff to precipitation, than Johnson for most months except March, which shows almost the same runoff ratio value in both catchments (see Fig. 2). Monthly runoff rates show the highest interbasin differences in the dry season (June-August). This is attributed to different infiltration mechanisms as well as to geographic characteristics such as slope, soil, and shape of the catchment. Due to different geology and soils, precipitation in Fanno is less infiltrated and rapidly reaches the river, while the more infiltrated precipitation in Johnson is evaporated in warm and dry climate conditions. In the wet season (November-March), continuing rainfall results in saturated soil condition that can behave like an impervious surface, so differences in the monthly runoff rate are smaller than those of the dry season. The coefficient of determination of daily streamflow between the two catchments also shows higher linear relations (above 0.74) for the wet season and lower relations (below 0.63) for the dry season (see Fig. 2).
Observed daily precipitation, maximum and minimum temperatures, and streamflow data are used for hydrologic modeling and downscaling of GCM simulations. The climate data are obtained from the National Oceanic and Atmospheric Administration Cooperative Observer Program  (NRCS, 1986), and land cover (PNWERC, 2002) are used.

Climate simulations and downscaling methods
Generally, the coupled atmosphere-ocean general circulation models (GCMs) are the best tools for projecting future climate in response to GHG emission forcing. GCMs have diverse horizontal and vertical grid resolutions, climate process description and approximation, parameterization of subgridscale phenomena, and initial condition (Randall et al., 2007). These different structures among GCMs cause the wide variations and biases in regional climate reproduction and projection (e.g., Im et al., 2011). Some GCMs fail to simulate regional inter-annual or decadal climate variability, which are important drivers of specific regional climate.
To estimate GCM performance in the Pacific Northwest, Mote and Salathé (2010) rank the 20 GCMs, implemented in IPCC AR4, based on 20th century bias, a global performance index (Achuta Rao and Sperber, 2006), and North Pacific variability of temperature, precipitation, and sea-level pressures (Mote and Salathé, 2010). The North Pacific variability represents the teleconnection effects of El Niño Southern Oscillation (ENSO) and Pacific Decadal Oscillation (PDO) and other large-scale climate processes over the Pacific Northwest (Hamlet et al., 2010). Based on the study of Mote and Salathé (2010), this study selects the three best GCMs, which are CNRM-CM3, ECHAM5/MPI-OM, and ECHO-G. Better GCM performance at simulating historical climate does not inevitably indicate a realistic projection under GHG forcing. However, if a GCM has poor performance for current important climate variability in the region, the derived regional changes for future should also be misleading (Prudhomme et al., 2002). No downscaling method can completely correct for the GCM's errors. Additionally, this approach provides some useful information such as weighted factor of GCM simulations (e.g., Tebaldi et al., 2005), and reducing ensemble numbers for future climate projection (e.g., Mote and Salathé, 2010).
This study uses two GHG emission scenarios, the A1B and B1 emission scenarios. Most global climate modeling groups generally employ A2, the A1B and B1 GHG emission scenarios (Randall et al., 2007) as high, medium and low emission scenarios for the 21st century, respectively. We focus on mid-century change for 2040-2069, in which period A2 and A1B show similar GHG emission forcing. Therefore, A1B and B1 emission scenarios can cover high and low GHG emission conditions. The climate simulation of three GCMs with two GHG emission scenarios are obtained from the World Climate Research Programme's (WCRP's) Coupled Model Intercomparison Project phase 3 (CMIP3) multi-model dataset (WCRP CMIP3, 2010).
To downscale three GCM simulations with two emission scenarios, we use a simple delta method, which has widely been used in climate change impact studies (e.g., Lettenmaier et al., 1999;Wilby and Harris, 2006;Loukas et al., 2007;Graham et al., 2007;Kay et al., 2009;Choi et al., 2009). This method first calculates monthly precipitation and temperature differences between the reference and future GCM simulations. Then, the obtained monthly differences between the two periods are applied to historical daily data for the reference period by adding monthly absolute differences for temperature and by multiplying percent differences for precipitation. This method can preserve the spatial and temporal variation of observation and remove the bias of GCM simulations. However, the delta method does not capture changes in precipitation and temperature variability from climate models and does not allow for more complex changes in daily extreme of precipitation and temperature (Hamlet et al., 2010). Therefore, changes in day-to-day variability of climate simulations are not taken into account in this study. This could lead to an underestimation of future flood frequency change.

Hydrologic model and parameter uncertainty
The PRMS model, Modular Modeling System (MMS) version developed by US Geological Survey (Leavesley and Stannard, 1996), is used in this study. This model simulates a water balance for each day and an energy balance for each half-day in each Hydrologic Response Unit (HRU), which is assumed to be homogeneous in its hydrologic response to given climate and land use conditions (Hay et al., 2009). A detailed description of the PRMS model structure is found in Leavesley et al. (2005). PRMS has seven parameters which are directly associated with land use change (see Table 1). Seasonal vegetation cover density (covden sum, covden win) and cover type (cov type) affect the amount of interception on HRUs. The seasonal vegetation cover density is determined by different leaf loss of cover types, such as grass, shrub, deciduous and coniferous trees (Viger and Leavesley, 2007, p. 99). Maximum values of interception storage for each cover type are considered by season and precipitation type (wrain intcp, srain intcp, and snow intcp). Ratio of impervious surface area on HRU (hru percent imperv) is a more important parameter in land use change impact on flood analysis, because it is highly sensitive to urbanization. High impervious surface area in this model induces less infiltration to soil and more overland flow to streams, potentially increasing peak flow volume.
PRMS is a physically-based hydrologic model, so some parameters can be obtained from physiographic characteristics and land surface features of the watershed using GIS layers, such as DEM, Land use, and Soil data . This study uses fixed parameters from GIS layers over time, except parameters related to land use. Snow effects are  Leavesley and Stannard (1996) for snow modeling in PRMS. We calibrate eight parameters that are associated with the timing and amount of runoff components (see Table 1). As in previous studies, streamflow simulation is most sensitive to these parameters (Bae et al., 2008b;Hay et al., 2009;Im et al., 2010;. LHS (McKay et al., 1979) is employed to sample the parameters from plausible ranges. LHS is an efficient sampling method that provides larger sample space with less computational effort comparable to those obtained from the conventional Monte Carlo simulation (Tang et al., 2007;Davey, 2008). LHS divides the feasible parameter space into equal intervals, so that at least one sample of each parameter set is selected randomly from each interval (Yang et al., 2010). To do an exhaustive search of behavioral parameters we decide to sample 20 000 parameters using LHS. These parameter sets are used to determine the closeness between daily simulated and observed streamflow for the period of 2000-2006 in both catchments. The Nash-Sutcliffe (1970) nondimensional model efficiency criterion (NSE) is used as a goodness of fit measure, with a value in excess of 0.6 indicating satisfactory fit between observed and simulated hydrographs (see Wilby, 2005;Choi and Beven, 2007). The NSE is generally more sensitive to high flow than low flow. Thus, an NSE score above 0.6 was considered appropriate for our flood frequency-focused study, since it mainly considers high flows. This approach can show the relative importance of parameter uncertainty in climate impact studies, although it cannot cover total equifinality of parameters (Beven, 2001). Therefore, the whole range of parameter uncertainty on flood frequency estimation is probably larger than what is presented in this study.

Natural variability
The natural variability of climate is the inherent internal fluctuation caused by combined effect of low-frequency (longer than 10 years) and high-frequency (shorter than 10 years) variability of nature (Wigley and Raper, 1990). The flood frequency analysis could be sensitive to the finite sampling within the natural climate variability (Kay et al., 2009). Therefore, it is essential that we analyze this sensitivity as an internal effect and anthropogenic climate change as external forcing. This will reveal the main source of uncertainty and indicate which source is a key controlling factor for future flood frequency change.
To estimate the effect of natural variability, this study applies simple replacement of climate time series using a moving block bootstrap method. The moving block bootstrap method (Künsch, 1989) is a resampling method with replacement to obtain a large number of samples (pseudo time series) for a time series, which have independent data structures such as precipitation (Ebtehaj et al., 2010). This Table 2. Description of the three land cover scenarios used in this study to simulated land cover projections within the Fanno and Johnson Creek catchments by 2050 (Source: Hulse et al., 2004;Franczyk and Chang, 2009 (Kay et al., 2009). For instance, the climate data of three months (December-February) in 1960 are randomly selected from any 3-month period between the water year 1960 and 1989. The selection of climate data with the same months is repeated 30 times until the years of new series are the same of original time series. This process allows the selection of data for a specific water year which could be repeated or may not be used at all. Flood frequency using 100 resampled climate series are compared to that obtained from original data. Also, the 100 resampled climate series are adjusted by the delta method described above to generate future climate conditions by the aforementioned three GCMs with two emission scenarios.

Flood frequency analysis -PeakFQ
To estimate the impacts of climate and land use changes on flood frequency, this study used typical statistical flood frequency analysis of maximum annual flood series using the PeakFQ program. PeakFQ provides estimates of instantaneous maximum annual peak-flows having diverse recurrence intervals such as 2, 5, 10, 25, 50, 100, 200, and 500 years as annual-exceedance probabilities of 0.50, 0.20, 0.10, 0.04, 0.02, 0.01, 0.005, and 0.002, respectively. Here, a 100 year flood describes a flood that is believed to have a probability of being equal or exceeding 0.01 in any one year (Raff et al., 2009). This program was developed based on the Bulletin 17B guidelines of the Interagency Advisory Committee on Water Data (IACWD, 1982), which is recommended for use by Federal agencies in the US. Bulletin 17B assumes that flood frequency can be described by a log-Pearson Type 3 (LP3) probability distribution (Griffis and Stedinger, 2007). Here, the LP3 distribution defines the probability that any single annual peak flow will exceed a specified streamflow. LP3 has three parameters: mean, standard deviation, and skew coefficient (Bobee and Ashkar, 1991). The skew coefficient is highly sensitive to the collected sample data of annual maximum floods, so that PeakFQ provides guidance on estimating the skew coefficient, such as the generalized skew from a digitized copy of the map in Bulletin 17B, the approach applied in this study.

Land use change scenarios
To consider possible future land use changes in both catchments, this study used two land-cover datasets developed by the Pacific Northwest Ecosystem Research Consortium (PN-WERC, 2002). The PNW-ERC provides three different land use scenarios for every 10 years of 2000-2050, namely, the conservation, the plan trend, and the development scenarios (see Table 2). These scenarios represent different future landscapes, based on projected human population growth patterns and potential development characteristics throughout the Willamette River basin (Hulse et al., 2004). As shown in Table 2, the conservation scenario assumes that greater emphasis on ecosystem protection and restoration will be implemented. The Plan Trend scenario assumes that current land use trends continue. The development scenario depicts greater expansion of urban growth boundaries (UGBs) with free rein to market forces across all components of the landscape, resulting in sprawl urban development. More detailed description of these scenarios is found in Hulse et al. (2004). This study used the conservation and the development scenarios as two extreme cases. A similar approach has been used in Franczyk and Chang (2009) and Praskievicz and Chang (2011).

Comparison of uncertainty sources
To identify the main source of uncertainty, we compare the maximum range of flood frequency change according to each uncertainty source . For instance, to determine the effect of GCM simulations (GCM structures), we first calculate the differences in flood frequency changes that are derived by different GCM simulations while holding the other data such as land use changes, emission scenarios, hydrologic model parameters, and natural variability constant. We then rank these differences and determine the maximum value at the top 5%. The same methodology is repeated to determine the maximum range for each uncertainty source.

Hydrologic model calibration
To calibrate PRMS model parameters, HRUs for the two catchments are delineated based on streamflow network, slope, aspect, and soil type. The geophysical parameters are extracted from DEM, land use, and soil GIS layers (see Table 1). The rest of the parameters (eight process parameters) are calibrated using Rosenbrock's (1960) automatic optimization method. The ratio of impervious surface area in HRU (hru percent imperv) is strongly related to land use change, as mentioned in Sect. 2.4. However, the land use layers of PNWERC do not provide the specific information of impervious surface area. They only describe some urbanrelated land use, such as residential, commercial, industrial, railroads, and roads. These land use categories contain both pervious and impervious surface areas. Therefore, if all urban land uses are assumed to be impervious surface areas, flood frequency would be overestimated. To determine the ratio of impervious surface area to urban land use, we develop an empirical relation between urban land use (%) and mean impervious surface area (%) (see Fig. 3) based on the data set of Waite et al. (2008). Waite et al. (2008) used different land use types, including mean impervious surface area, for 28 catchments in Oregon and Washington to estimate the effect of urbanization on steam ecosystems. As shown in Fig. 3, the estimated regression equation shows a good fit between urban land use and mean impervious surface area (R 2 = 0.99). The regression coefficients are used to estimate percent impervious surface areas in each HRU (hru percent imperv) in PRMS modeling for these two urban catchments.

Projected future climate change and land use change
Changes in monthly precipitation show different patterns by GCMs and GHG emission scenarios, but the changes are similar in the two catchments (see Fig. 4). The CNRM-CM3 and the ECHAM5/MPI-OM simulations project slight  Table 1). increases in winter (December, January, and February) precipitation, while predicting drier summers (June, July, August, and September) as indicated by previous studies (e.g., Mote et al., 2003;Graves and Chang, 2007;. In the study catchments, winter precipitation is closely related to flood events. Therefore, rising water tables resulting from an increase of winter precipitation and soil moisture content are likely to lead to more frequent flooding in this region. However, the ECHO-G projects a slight decrease in winter precipitation. These different precipitation projections contribute to uncertainty in flood frequency analysis. Climate change projection for monthly temperatures ranges from +0.3 • C increase in February (CNRM-CM3, B1) to +6.1 • C in August (ECHO-G, A1B) for the 2050s (not shown). Figure 5 shows changes in land use categories of three different land use data sets -reference land use in 2001, the conservation and the development land uses for the 2050s. The two catchments are projected to have different paths of future growth, as reflected in changes in each land use category. In the Fanno creek catchment, absolute changes in land use categories are small because it is already highly developed (85% in 2001). Hence, the Johnson creek catchment shows considerable differences in each land use among the three scenarios. Urban land use shows a 17% increase under the development (sprawl development) scenario and an 11% increase under the conservation (compact development) scenario because of population growth, construction of building and roads, and urban development in agricultural land use (Hulse et al., 2004). Agricultural land use in both future scenarios decreases by approximately 17%. Grass-land and forest land uses are higher under the conservation scenario than under the development scenario.  Figure 6 shows the range of flood frequency at the reference and future climate change conditions, excluding land use change effects. The reference period only considers the natural variability impact. Hence, the two future periods represent impacts of climate change on flood frequency caused by the combined conditions of climate change and natural variability. The effect of climate change is much more dominant in both catchments as compared with natural variability (taller box and whisker). The t-test results show that the flood frequency of all return periods significantly changes due to climate change at the 95% confidence level (see Table 3). The GHG emission scenarios are only significantly different for 2-year flood frequency. The climate change impact on flood frequency between both catchments is similar. This is attributed to the fact that the catchments are located in same climate region in the Willamette Valley and analyses are made using data derived from coarse scale GCM simulations. In a contrasting case study, Kay et al. (2009) show different responses between two distant catchments in UK using regional climate model (RCM) simulations. They show that one catchment is highly dominated by natural variability, while the other catchment is strongly affected by climate change. Hulme et al. (1999) explain that if a region is more dominated by natural variability than by climate change, adaptation management that takes into account natural variability may be sufficient to withstand climate change. Our results show that future flood management in the Fanno and Johnson creek catchments should consider climate change impact as well as historical natural climate variability. As shown in Fig. 7, the natural variability impact is much greater than future land use change impact. The variation in flood frequency caused by land use change is similar to that due to natural variability in both catchments. However, under the development scenario, short-term floods ( confidence level (see Table 3). This might indicate that notable land use change in less developed catchment could lead to significantly more frequent bankfull flooding, although the natural variability effect is more pronounced for larger flood events. The median values of flood frequency under the development condition are slightly higher than those of the conservation scenario. Also, shorter term floods increase more than longer term floods.

Projected flood frequency
For the combined impact of climate and land use changes, flood frequency at the six different return periods increases slightly, though each change had high variations (Fig. 8). The range of flood frequency change gradually increases from shorter term floods to longer term floods. The variations under the A1B scenario are larger than those under the B1 scenario in both catchments. Since variation is high, an interpretation of the flood frequency impact of each scenario solely based on Fig. 8 is difficult. Accordingly, we calculated ensemble mean value of flood frequency change for each scenario. Figure 9 shows the ensemble mean of relative changes of flood frequency under two GHG emission, two land cover change, and the combined scenarios (four) that are calculated from the reference flood frequency. The A1B scenario shows the biggest change among the separate emission and land cover scenarios in both catchments. In the Fanno creek catchment, ensemble results of all 8 scenarios show higher changes than those caused by natural variability. However, in Johnson creek, the natural variability impact becomes more significant than the B1 and land cover change scenarios for short-term flood frequency of less than 25 year floods. In all cases, the combined impacts on flood frequency are higher than those of natural variability in both catchments. Of the combined land use and climate scenarios, the A1B with development scenario induces the highest increase in flood frequency, and the B1 with conservation scenario induces the lowest increase in flood frequency. The shorter term flood frequencies are more sensitive to the combined scenarios than longer term ones (see Table 4). Further, the difference between A1B with development scenario and B1 with conservation scenario is greater in Johnson than in Fanno (see % difference between the two scenarios in Fig. 9). For the long term extremes, the Johnson creek shows significant difference between the A1B with development scenario (6.6% difference) and the B1 with conservation scenario (3.4% difference) (see Table 4).
This result indicates that, under expected climate change conditions, an adaptive urban planning based on the conservation scenario could be more effective in less developed Johnson catchment than in the already developed Fanno. Also, this result demonstrates that the combined effect does not linearly increase catchment flood frequency. For example, 2 years floods in Fanno are increased by 12.4% by the A1B scenario alone, and by 9.7% by the development scenario alone; however, they are increased by 14.8% by the combination of the A1B and development scenarios. This could be attributed to nonlinear hydrologic responses under different climate and land use conditions. Additionally, it implies that if we want to obtain more realistic future    projections of urban flood risk, we need to develop possible climate change scenarios as well as land use change scenarios. Figure 10 shows the relative variation (uncertainty) in flood frequency change projections under the combined impact of climate and land use change. Uncertainty due to land use change is the smallest in this study, except for the occurrence of 2 year floods at Johnson creek, although the Johnson's range is larger than the Fanno's. This could indicate that longer term floods could be less affected by land use change than climate change. However, this result also suggests that if land use at a catchment scale changes abruptly, the land use change will become a more significant uncertainty source than climate change for short term floods. Emission scenario uncertainty also shows a relatively smaller range than those of the other sources. The uncertainty from hydrologic parameters is more significant at Fanno than Johnson, but it is smaller than uncertainty due to GCM and natural variability. GCM uncertainty strongly affects shorter term 2 and 5 year floods, while longer term 25, 50, and 100 year floods are more controlled by natural variability. This demonstrates that both uncertainty sources, GCMs and natural variability, are significant factors in urban flood frequency analysis.

Caveats of this study
This research deals with uncertainty in future flood frequency analysis in two distinct urban areas. We consider three uncertainty sources; climate projection (GCM structure, future The vertical ranges show the 95% confidence interval. 4 GHG emission scenario, and natural variability), urban development (different future land use planning scenarios), and hydrologic modeling (hydrologic model parameters). Our results contribute to an understanding of the combined effects of hydro-climatic modeling and urban development effects on urban flood analysis. While we identify the relative magnitude of uncertainties arising from the sources mentioned above, there are remaining uncertainty sources, such as GCM initial condition, downscaling method, and hydrologic model structure, which are not investigated in the current study. Therefore, our results should be cautiously interpreted along with other potential sources of uncertainties. We carefully select the three best GCMs, but these GCMs do not necessarily project future climate accurately. Furthermore, three GCMs are insufficient to cover the full range of GCM structure uncertainty. However, our results show the uncertainty caused by GCMs is higher than that due to other sources. This is consistent with the findings of previous studies (e.g., Wilby and Harris, 2006;Kay et al., 2009). Therefore, the end-to-end effect of GCM uncertainty on flood frequency projection could be larger than that presented in this study, however, the relative magnitudes of the GCM structural uncertainty might not vary. The uncertainties due to future GHG emissions are not fully considered as proposed in the IPCC storyline (IPCC, 2000).
Our results are also affected by the simple delta method for downscaling GCMs because this approach cannot consider changes in interannual or day-to-day variability of climate simulations (Im et al., 2010;Prudhomme and Davies, 2009).
Additionally, using a different NSE threshold value could have resulted in a wider or narrower parameter uncertainty range, although it would still not be significant compared to other uncertainty sources. However, more sophisticated methods and approaches in quantifying the parameter uncertainty relying on Sequential Monte Carlo (SMC) using ensemble filtering (Moradkhani et al., 2005a,b;Moradkhani and Sorooshian, 2008;Leisenring and Moradkhani, 2010;DeChant and Moradkhani, 2011;Montzka et al., 2011), Markov Chain Monte Carlo (MCMC) (e.g., Smith and Marshall, 2008;Vrugt et al., 2009) and Moving Block Bootstrap Sampling (MBBS) (Ebtehaj et al., 2010) can be employed. Furthermore, we did not include uncertainties associated with climate data downscaling (Fowler et al., 2007;Im et al., 2010;Najafi et al., 2010a), and hydrologic model structure (Clark et al., 2008;Jiang et al., 2007;Bae et al., 2011;Najafi et al., 2010b). Recently, Najafi et al. (2010b) used Bayesian Model Averaging to quantify and minimize the uncertainty associated with hydrologic model structure and selection in the context of hydrologic climate change impact studies.
Urban climate is controlled not only by global and regional natural climate systems, but also by local urbanization effects, such as the urban heat island, the urban canopy layer, and varying aerosol composition (Ntelekos et al., 2010). Urbanization could significantly affect the precipitation climatology relating to flood events (Shepherd, 2005). Ntelekos et al. (2008) demonstrates that rainfall accumulations of 30% of the total extreme events are attributed to urbanization impacts in the Baltimore metropolitan area, Washington DC. Therefore, the interaction between global climate change and urban climatology is another important uncertainty source in urban climate impact studies.
In changing climate conditions, an assumption of stationarity in flood frequency analysis may not be valid (Milly et al., 2008;Smith et al., 2005). This study uses the PeakFQ based on the Bulletin 17B that assumes the constant distribution of flood events regardless of climate change. Some previous studies illustrate that a traditional approach to flood frequency estimation could not rely on stationarity assumptions (Raff et al., 2009;Sivapalan and Samuel, 2009). Now, a robust methodology for incorporating projected climate information into flood frequency analysis is needed.

Conclusions
This study examines the potential changes of flood frequency and the associated uncertainties in the two catchments exhibiting different levels of urbanization. Here, the important conclusions are summarized.
1. For the combined scenarios, GCM uncertainty highly affects shorter term extremes, while longer term extremes are more controlled by natural variability. Hence, the uncertainties due to future GHG emission scenarios and land use change scenarios are less important than natural variability. Also, hydrologic model parameter has less impact than natural variability and GCM structure in our uncertainty analysis.
2. The combined impacts of land use change and climate change scenarios induce significant changes in the shorter term extremes in both catchments. Flood frequency change demonstrates the highest increase under the A1B with development scenario and the lowest increase under the B1 with conservation scenario.
3. In the 2050s period, flood frequency is projected to slightly increase in both catchments, although there are substantial uncertainties. Changes in flood frequency are more sensitive to climate change (A1B scenario) than land use change. Land use change impact is only significant in the less developed Johnson catchment, which is projected to be more urbanized in the 2050s.