Flood Risk Assessment in Urban Areas of Southern Taiwan

A flood risk assessment of urban areas in Kaohsiung city along the Dianbao River was performed based on flood hazards and social vulnerability. In terms of hazard analysis, a rainfallrunoff model (HEC-HMS) was adopted to simulate discharges in the watershed, and the simulated discharges were utilized as inputs for the inundation model (FLO-2D). Comparisons between the observed and simulated discharges at the Wulilin Bridge flow station during Typhoon Kongrey (2013) and Typhoon Megi (2016) were used for the HEC-HMS model calibration and validation, respectively. The observed water levels at the Changrun Bridge station during Typhoon Kongrey and Typhoon Megi were utilized for the FLO-2D model calibration and validation, respectively. The results indicated that the simulated discharges and water levels reasonably reproduced the observations. The validated model was then applied to predict the inundation depths and extents under 50-, 100-, and 200-year rainfall return periods to form hazard maps. For social vulnerability, the fuzzy Delphi method and the analytic hierarchy process were employed to select the main factors affecting social vulnerability and to yield the weight of each social vulnerability factor. Subsequently, a social vulnerability map was built. A risk map was developed that compiled both flood hazards and social vulnerability levels. Based on the risk map, flood mitigation strategies with structural and nonstructural measures were proposed for consideration by decision-makers.


Introduction
In recent years, the economy has developed rapidly in various parts of the world, and countries have used many fuels, such as oil and coal, that cause holes in the ozone layer and produce a large amount of greenhouse gases such as carbon dioxide and methane, resulting in a continuous rise in temperature. Droughts, floods, and other extreme weather events occur continuously, and extreme weather causes natural disasters that continue to threaten human lives and cause severe property losses [1][2][3][4][5]. According to statistics from the United Nations Disaster Database (www.emdat.be, (accessed on 30 January 2021)), since 1950, the number of major natural disasters has increased year by year, and after 1990, it has exceeded 200 annually. The total number of disasters that occurred in 2018 was 282. In the seven continents in the world, the number of disasters that occurred in Asia in 2018 was 129, accounting for 46% of the total. Among the types of natural disasters that occurred in 2018, flood-related disasters accounted for 74% of all disasters; floods made up 39% of the total number of disasters, storms made up 30%, and slope disasters made up 5%. The rest of the disaster types account for approximately 2% to 9% each. Therefore, quantifying flood risk is an important way to develop adaptation and management strategies for decision-making [6,7]. In general, flood risk can be defined as the product of flood hazard and vulnerability [8][9][10].
Based on a literature review, several approaches have been utilized to perform flood risk assessments. These methods include historical disaster flood risk assessments [11][12][13][14][15][16][17], index systems for flood risk assessments [18][19][20][21][22][23], scenario simulation approaches [24][25][26][27][28][29][30], and remote sensing and geographic information systems [31,32]. Different approaches display corresponding advantages and disadvantages [26]. Cai et al. [26] employed a one-dimensional pipe network and a two-dimensional surface-coupled model and a multiindex fuzzy comprehensive evaluation model to build flood disaster risk. Wang et al. [23] utilized a fuzzy synthetic evaluation method based on combined weight to generate flood risk maps in the Beijing-Tianjin-Hebi urban area. Geng et al. [33] proposed a coupled oneand two-dimensional hydrodynamic model to simulate inundation as a hazard indicator, and social and economic characteristics were analyzed as vulnerability indicators. The analytic hierarchy process (AHP) was used to compute the weights of the assessment indices. Afterwards, a fuzzy comprehensive evaluation approach was adopted to assess urban flood risk. However, these studies did not clearly describe that how the social vulnerability was selected for further analysis.
Taiwan is located on the western edge of the Pacific Ocean, surrounded by the sea; it is an island country. The river slopes are steep and winding. It is often exposed to extreme rains and typhoons during the rainy season and in summer and autumn. In the wet season, southwesterly winds and typhoons bring frequent heavy rains, resulting in failures of regional drainage systems in the middle and lower reaches of rivers; these failures cause levees to break or overflow, resulting in flooding and threatening the lives of the people in the middle and lower reaches of the rivers. The property and economic losses associated with these floods are serious. In the future, rainfall patterns will change as a result of climate change. There may be more rainfall during the wet season, and flooding situations in urban regions and low-lying areas will become more serious. Therefore, the government must take preventive measures with flood control projects such as constructing urban drainage systems, riverbanks, and rainwater sewers as early as possible to ensure the safety of people and property.
In addition to establishing flood warning systems and improving existing drainage systems, governments are now gradually moving toward flood prevention and disaster mitigation planning; these goals are realized through the formulation of laws and regulations, the implementation of disaster preventative measures, and mitigation education and training aimed at achieving effective risk reduction. Flood risk assessments are composed of hazards and vulnerability. Among these, hazards are associated with occurrences of natural disasters that humans cannot eliminate or change through existing technology. To effectively reduce risk, we must decrease the level of vulnerability.
The objective of the present study is to utilize the hydrological model HEC-HMS and the inundation model FLO-2D to explore the areas that are prone to flooding in the urban region of the Dianbao River watershed in southern Taiwan and draw a flood potential map using the degree of hazard. The fuzzy Delphi method and the hierarchical analysis approach were utilized to analyze and calculate the degree of vulnerability. By combining the degree of hazard and the degree of vulnerability, an urban area flood disaster risk map was produced. Based on the risk map, adaptation strategies were proposed for disaster prevention and mitigation planning.
The hypotheses of this study include that (1) the minimum unit for flood risk assessment is village. It means that flood risk level in a village is same; (2) no engineering works is executed to reduce the flood hazard; and (3) the social vulnerability for different rainfall return periods is also same.

Data Collection
This study collected rainfall data at eight rain gauge stations, flow data at the Wulilin Bridge flow station, and water level data at the Changruu water level station (see Figure 1a). The collection period of the rainfall, flow, and water level data was from August 2012 to March 2018. of this pond is approximately 425,000 m 3 , while the area of the flood detention pond in Area B is 42 hectares and its flood detention volume is approximately 1.05 million m 3 . In this study, the data of rainwater sewers and flood detention ponds in the middle and lower reaches of the Dianbao River were taken from the Taiwan Experimental Watershed Information Platform (Figure 1d).
Apart from the data collection used for the rainfall-runoff model and inundation model, data related to social vulnerability were also collected for further analysis.  The land use data were obtained from the Taiwan Experimental Watershed Information Platform established by the Taiwan Typhoon Flood Research Institute and the Water Resources Agency. The land use types can be classified as agriculture, forest, traffic, hydraulic, construction, public, recreation, mineral, or others ( Figure 1b).
For the elevation data of the Dianbao River watershed, 20-m DTM data released from the information disclosure platform of the Ministry of the Interior were adopted (Figure 1c). A model was utilized to generate grids using DTM data to obtain the numerical elevation of the terrain, which was further used to determine the accuracy of the model simulations. The cross-sectional river data were collected from the Taiwan Experimental Watershed Information Platform.
In urban areas, cities are frequently flooded due to heavy rainfall events with short durations and the rapid collection of water resulting from urban development. Therefore, the construction of rainwater sewer systems and flood detention ponds would help urban areas drain excessive rainwater quickly. In the Dianbao River watershed, there are flood detention ponds in zones A and B, located on both banks of the Daliao drainage area. The area of the flood detention pond in zone A is 17 hectares, and the flood detention volume of this pond is approximately 425,000 m 3 , while the area of the flood detention pond in Area B is 42 hectares and its flood detention volume is approximately 1.05 million m 3 . In this study, the data of rainwater sewers and flood detention ponds in the middle and lower reaches of the Dianbao River were taken from the Taiwan Experimental Watershed Information Platform (Figure 1d).
Apart from the data collection used for the rainfall-runoff model and inundation model, data related to social vulnerability were also collected for further analysis.

Study Area
The boundaries of the Dianbao River watershed in Kaohsiung city border the Agongdian River watershed in the north, the Houjin River drainage basin in the south, the hilly area of the Central Mountain Range in the east, and Mituo District in the west. The basin originates from the top of Wushan Mountain in Yanchao District, Kaohsiung city. The elevation is approximately 320 m, and the river flows westward through Yanchao District, Dashe District, Nanzi District, Qiaotou District and Gangshan District and finally into the Taiwan Strait on the south side of the Oziliao fishing port in Ziguan District [34]. Most of the Dianbao River watershed consists of plains. The main roads passing through the area are high-speed highways, longitudinal railways, other highways and roads, such as county roads. The road network is quite dense, and the traffic is convenient for users. This study uses the Taiwan High Speed Rail as a boundary to simulate rainfall-runoff with the HEC-HMS model in the upstream watershed area in the east and to simulate the flooding area with the FLO-2D in the west (Figure 1a).
The average annual rainfall in the watershed of the Dianbao River is approximately 1825 mm [35]. The wet season is mainly from June to August, and the dry season is from October of a given year to March of the following year.

Methodological Outline
The procedure followed in this research includes two parts: hazards and vulnerability. Before starting the research, fundamental data, including terrain and land classification, urban drainage networks, basin geography, meteorological and hydrological observations, population, social economy, rescue equipment, and special agency data, were gathered. In terms of hazards, the research steps were to establish the hydrological and inundation models, perform model calibration and validation, predict the flooding extents and depths under different rainfall return periods, and finally complete the hazard map. In terms of social vulnerability, the research steps were to select social vulnerability factors, design questionnaires, apply the fuzzy Delphi method to determine the main factors affecting social vulnerability, use the AHP to determine the weights of the main social vulnerability factors, calculate the social vulnerability index (SVI), and finally establish a social vulnerability map. According to these two maps, the risk map was built. The research outline is illustrated in Figure 2.

Study Area
The boundaries of the Dianbao River watershed in Kaohsiung city border the Agongdian River watershed in the north, the Houjin River drainage basin in the south, the hilly area of the Central Mountain Range in the east, and Mituo District in the west. The basin originates from the top of Wushan Mountain in Yanchao District, Kaohsiung city. The elevation is approximately 320 m, and the river flows westward through Yanchao District, Dashe District, Nanzi District, Qiaotou District and Gangshan District and finally into the Taiwan Strait on the south side of the Oziliao fishing port in Ziguan District [34]. Most of the Dianbao River watershed consists of plains. The main roads passing through the area are high-speed highways, longitudinal railways, other highways and roads, such as county roads. The road network is quite dense, and the traffic is convenient for users. This study uses the Taiwan High Speed Rail as a boundary to simulate rainfall-runoff with the HEC-HMS model in the upstream watershed area in the east and to simulate the flooding area with the FLO-2D in the west (Figure 1a).
The average annual rainfall in the watershed of the Dianbao River is approximately 1825 mm [35]. The wet season is mainly from June to August, and the dry season is from October of a given year to March of the following year.

Methodological Outline
The procedure followed in this research includes two parts: hazards and vulnerability. Before starting the research, fundamental data, including terrain and land classification, urban drainage networks, basin geography, meteorological and hydrological observations, population, social economy, rescue equipment, and special agency data, were gathered. In terms of hazards, the research steps were to establish the hydrological and inundation models, perform model calibration and validation, predict the flooding extents and depths under different rainfall return periods, and finally complete the hazard map. In terms of social vulnerability, the research steps were to select social vulnerability factors, design questionnaires, apply the fuzzy Delphi method to determine the main factors affecting social vulnerability, use the AHP to determine the weights of the main social vulnerability factors, calculate the social vulnerability index (SVI), and finally establish a social vulnerability map. According to these two maps, the risk map was built. The research outline is illustrated in Figure 2.

Hazard Assessment with Hydrological and Hydrodynamic Models
In recent years, the impacts of climate change have caused extreme rainfall events, which have resulted in urban areas being unable to discharge internal waters, leading to flooding. This study uses the hydrological model HEC-HMS and the inundation model FLO-2D to simulate an inundation map of urban areas, serving as the hazard map.

Hydrological Model (HEC-HMS)
The rainfall-runoff model used in this study, HEC-HMS software, was developed by the Hydrological Engineering Center of US Army Corps Engineers. The HEC-HMS model, which was modified from HEC-1, is a new version used for hydrological simulations [36]. In this study, HEC-HMS version 4.3 was employed to calculate the rainfall-runoff discharge in the watershed [36]. The model was designed to be applied to different geographic areas, such as in studies of river basin flood hydrology and large and small watershed/urban runoff [37][38][39][40][41][42][43]. The HEC-HMS includes three modules: the basin model, meteorological model, and control specification. In the basin model, the US Soil Conservation Service (SCS) curve number (CN) was used to calculate precipitation loss. The input data include the initial loss, CN, and impermeability. The HEC-GeoHMS toolbox can be adopted to calculate the average CN for each subwatershed. The SCS synthetic unit hydrograph method was applied to convert rainfall to runoff discharge. Two parameters, the time of concentration (T C ) and storage coefficient (R), were used in the Clark unit hydrograph method. The time of concentration (T C ) and the regional value (RV), which is related to the time of concentration and the storage coefficient, are expressed as follows [44,45]: where L denotes the length of the main channel (m) and S expresses the slope of the watershed (%). The time lag (T lag ) then can be determined as follows: The recession base flow method, which uses an exponentially decayed algorithm, was employed to separate the base flow. In a rainfall event, the recession base flow strongly affects the volume and timing of the base flow. In this study, the recession constant and ratio to the peak were set in the ranges of 0.2-0.35 and 0.01-0.1, respectively.
The rainfall method was used to establish a meteorological model in this study. In the module, the time-series rainfall data of each subcatchment and the total rainfall of each subcatchment are given. The main parameters to be set in the control specification are the start time and end time of the simulation and the time step. In this study, the computational time step was set to one hour.

Inundation Model (FLO-2D)
FLO-2D model is generally adopted for studies on flood routing. The model is currently certified by the Federal Emergency Management Agency. It is mainly used in urban flood simulations, flood plain management, engineering risk designs, and debris flow simulations. FLO-2D is a combination of a one-dimensional river channel unsteady flow model and a two-dimensional horizontal overland flow model [46]. It is used to simulate the depth and velocity of a two-dimensional fluid flow and to estimate the reasonable disaster range of a fluid. The detailed governing equations can be found in Dimitriadis et al. [47] and Hu et al. [48].
The FLO-2D model divides the DEM (digital elevation model) into fixed-size square grids. Each grid can be input with a surface elevation value, Manning roughness coefficient, Sustainability 2021, 13, 3180 6 of 22 infiltration coefficient, area reduction factor, outflow, inflow, and other hydrological and geological data to simulate the real flood flow in the study area. This model starts with basic two-dimensional surface overland flow, and many data parameters can be added; users can freely choose to turn on or off the additional functions, including the presence of river channels, dikes, roads, streets, urban area reduction, rainfall, debris flow, infiltration, and rainwater sewers.

Establishment of HEC-HMS and FLO-2D Models
To save computational time, the simulation domain of the FLO-2D model is smaller than that of the HEC-HMS model; therefore, connecting these two modeling domains is necessary. Figure 3a illustrates the connection of the HEC-HMS and FLO-2D models. It shows that the whole domain of the Dianbao River watershed is used for HEC-HMS, while the urban area is adopted for FLO-2D.
Using geomorphological analysis with the input function of the HEC-HMS model and the output function of HEC-GeoHMS can achieve the purpose of automatically building watersheds. The Dianbao River watershed was built for the HEC-HMS model as shown in Figure 3b and includes 7 subcatchment units, 16 channel units, and 16 confluence point units.
To simulate the inundation extent and depth of urban areas, a numerical grid was established in the FLO-2D model. The DTM data were imported into the FLO-2D model, and using the grid generation module, a grid with a size of 20 m × 20 m in the horizontal plane was generated. The new version of FLO-2D has been enhanced and can be combined with the EPA-SWMM model. Therefore, the grid system calculations can be used to simulate the drainage of rainwater sewer nodes and the overflow of full pipes. Figure 3c displays the grid generation of the FLO-2D model, which contains a total of 142,608 grids. The FLO-2D model divides the DEM (digital elevation model) into fixed-size square grids. Each grid can be input with a surface elevation value, Manning roughness coefficient, infiltration coefficient, area reduction factor, outflow, inflow, and other hydrological and geological data to simulate the real flood flow in the study area. This model starts with basic two-dimensional surface overland flow, and many data parameters can be added; users can freely choose to turn on or off the additional functions, including the presence of river channels, dikes, roads, streets, urban area reduction, rainfall, debris flow, infiltration, and rainwater sewers.

Establishment of HEC-HMS and FLO-2D Models
To save computational time, the simulation domain of the FLO-2D model is smaller than that of the HEC-HMS model; therefore, connecting these two modeling domains is necessary. Figure 3a illustrates the connection of the HEC-HMS and FLO-2D models. It shows that the whole domain of the Dianbao River watershed is used for HEC-HMS, while the urban area is adopted for FLO-2D.
Using geomorphological analysis with the input function of the HEC-HMS model and the output function of HEC-GeoHMS can achieve the purpose of automatically building watersheds. The Dianbao River watershed was built for the HEC-HMS model as shown in Figure 3b and includes 7 subcatchment units, 16 channel units, and 16 confluence point units.
To simulate the inundation extent and depth of urban areas, a numerical grid was established in the FLO-2D model. The DTM data were imported into the FLO-2D model, and using the grid generation module, a grid with a size of 20 m × 20 m in the horizontal plane was generated. The new version of FLO-2D has been enhanced and can be combined with the EPA-SWMM model. Therefore, the grid system calculations can be used to simulate the drainage of rainwater sewer nodes and the overflow of full pipes. Figure 3c displays the grid generation of the FLO-2D model, which contains a total of 142,608 grids.

Fuzzy Delphi (FD)
The Delphi method proposed by Dalkey and Helmer [49] is a procedural method used to express expert opinions in a systematic manner. However, the Delphi method has some shortcomings; it is time-consuming, and the cost of collecting expert opinions can be high. The so-called consistency of expert opinions only falls within a certain range of expert opinions but does not determine the point of impact. Therefore, this law implies ambiguity, so the final results cannot truly reflect the opinions of experts. To resolve these shortcomings, the fuzzy Delphi method was developed. Murray et al. [50] first combined the fuzzy set with Delphi, and then Ishikawa [51] proposed the max-min method and the fuzzy integration method to integrate the opinions of experts into fuzzy numbers, which is called the fuzzy Delphi (FD) method. In the first stage of the expert questionnaire, the method proposed by Jeng [52] was adopted to select important social vulnerability factors.
The expert questionnaire is mainly filled out based on the judgment of each expert's professional literacy and the importance of each evaluation factor. The evaluation scale is divided into levels from 0-10. The higher the evaluation score is, the more representative the factor is of the research trend. An option column was added at the end of the questionnaire. If the experts felt that there were deficiencies, their opinion can be supplemented, making the questionnaire more complete and in line with the meaning of the fuzzy Delphi questionnaire. The consensus value can be yielded by triangular fuzzy number method [50]. The results using FD to obtain important social vulnerability factors can be found in Section 4.4.1.

Analytic Hierarchy Process (AHP)
Saaty [53] proposed the analytic hierarchical process (AHP) method. Through the establishment of a hierarchical structure with mutual influence, the AHP can be used to solve multiple evaluation criteria and decision-making problems associated with uncertainty; through quantitative judgments and comprehensive evaluation, the AHP can provide decision-makers with sufficient information to choose appropriate options and reduce the risk of making mistakes. The purpose of the development of the hierarchical analysis method is to systematize complex issues, provide hierarchical decomposition from different levels, organize the levels in a quantitative manner, and comprehensively evaluate them to provide decision-makers with information with which they can choose appropriate solutions.
The AHP builds complex multi-objective decision-making problems into hierarchical systems with dendritic structures; each layer is composed of different elements. The

Fuzzy Delphi (FD)
The Delphi method proposed by Dalkey and Helmer [49] is a procedural method used to express expert opinions in a systematic manner. However, the Delphi method has some shortcomings; it is time-consuming, and the cost of collecting expert opinions can be high. The so-called consistency of expert opinions only falls within a certain range of expert opinions but does not determine the point of impact. Therefore, this law implies ambiguity, so the final results cannot truly reflect the opinions of experts. To resolve these shortcomings, the fuzzy Delphi method was developed. Murray et al. [50] first combined the fuzzy set with Delphi, and then Ishikawa [51] proposed the max-min method and the fuzzy integration method to integrate the opinions of experts into fuzzy numbers, which is called the fuzzy Delphi (FD) method. In the first stage of the expert questionnaire, the method proposed by Jeng [52] was adopted to select important social vulnerability factors.
The expert questionnaire is mainly filled out based on the judgment of each expert's professional literacy and the importance of each evaluation factor. The evaluation scale is divided into levels from 0-10. The higher the evaluation score is, the more representative the factor is of the research trend. An option column was added at the end of the questionnaire. If the experts felt that there were deficiencies, their opinion can be supplemented, making the questionnaire more complete and in line with the meaning of the fuzzy Delphi questionnaire. The consensus value can be yielded by triangular fuzzy number method [50]. The results using FD to obtain important social vulnerability factors can be found in Section 4.4.1.

Analytic Hierarchy Process (AHP)
Saaty [53] proposed the analytic hierarchical process (AHP) method. Through the establishment of a hierarchical structure with mutual influence, the AHP can be used to solve multiple evaluation criteria and decision-making problems associated with uncertainty; through quantitative judgments and comprehensive evaluation, the AHP can provide decision-makers with sufficient information to choose appropriate options and reduce the risk of making mistakes. The purpose of the development of the hierarchical analysis method is to systematize complex issues, provide hierarchical decomposition from different levels, organize the levels in a quantitative manner, and comprehensively evaluate them to provide decision-makers with information with which they can choose appropriate solutions.
The AHP builds complex multi-objective decision-making problems into hierarchical systems with dendritic structures; each layer is composed of different elements. The complex problems are changed from high-level problems to low-level problems as they are gradually decomposed, and each level in the framework only affects one other level and, at the same time, is only affected by one other level. Before establishing the AHP framework, one thing must be confirmed: each factor must be independent of the other factors. The hierarchical structure method can be used to confirm a hierarchical relationship, but there is no specific construction method in the actual construction process. The AHP uses the characteristic vector method to calculate the weights of the elements and obtains the priority value of each scheme. The larger the value is, the higher the priority of the adopted scheme is.
Therefore, this study adopts the AHP method to design expert questionnaires, which are completed by experts in related fields, and the weights of each vulnerability factor are determined based on the results of the questionnaire. The results using AHP to obtain the weight of social vulnerability factors can be found in Section 4.4.2.

Social Vulnerability Index (SVI) and Normalization
This study uses the FD method to select and screen out the vulnerability factors and then uses the AHP method to determine the weight of each vulnerability factor. In terms of calculating vulnerability, this research refers to the social vulnerability index proposed by Shaw et al. [54] to identify the vulnerability of a certain area. Then, the vulnerability data are standardized based on this index. Since the values of the variables are different and the units are also different, to compare different variables with each other, it is necessary to standardize all the statistical data obtained before the vulnerability index is calculated. The index is calculated in an average manner. The standardized formula is expressed as follows: where x denotes the statistics of different variables, M represents the mean value (= 1 n expresses the total number of variables, SD indicates the standard deviation, and α is the correction factor of each vulnerability factor (=+1 or −1). The correction factor (α = −1) decreases the level of social vulnerability, representing a positive effect on social vulnerability, whereas the corrector factor (α = +1) increases the level of social vulnerability, denoting a negative effect on social vulnerability. The social vulnerability index (SVI) can be expressed as follows: where k is the weight value, Z denotes the standardized value, and N indicates the total number of vulnerability factors. The next step is to normalize the SVI value and convert it to a value between 0 and 1 (=I SVI ). When the I SVI value is closer to 1, the vulnerability is higher; in contrast, when the I SVI value is closer to 0, the vulnerability is lower. The normalization formula can be represented as follows: where SVI m represents the minimum SVI value and D denotes the full range of SVI values, which is the maximum SVI value minus the minimum SVI value.
After calculating the I SVI value and dividing the social vulnerability into 5 grades, including extremely low, low, medium, high, and extremely high (see Table 1) based on every 0.2 increment [55], ArcGIS software was utilized to draw the social vulnerability map.

Hazard Assessment
Urban flooding disasters are mainly caused by heavy rainfall. When rainfall exceeds the rainfall volume that the drainage system of an area can load, flooding disasters occur, causing casualties and property losses and indirectly leading to the suspension of commercial activities and classes in schools. In this study, the simulated inundation depth is used as the hazard assessment index, and the inundation depth is divided into 5 grades: below 0.3 m, 0.3 m to 0.5 m, 0.5 m to 1.0 m, 1.0 m to 3.0 m and above 3.0 m (see Table 1). When the inundation depth is 0.3 m, the water depth is approximately the height of a child's knees, so it is difficult for children to escape. When the inundation depth is 0.5 m, the water depth is approximately the height of an adult's knees, which causes adults to have difficulty walking. When the inundation depth is 1.0 m, the water depth is approximately the height of the waist of an adult, so it indeed causes difficulties for adults walking. When the inundation depth is 3.0 m, the water depth is as high as a building is tall, which threatens human lives.

Risk Assessment
Different definitions of risk have been documented by various researchers [56,57]. Generally, risk can be defined as the product of hazard (H) and social vulnerability (V) in Taiwan [58]. The simple equation to describe the risk is expressed as follows: Because hazards and social vulnerability are each graded in five levels, the approach of a semiquantitative risk assessment can be defined as the risk (Water Resources Agency, 2015). The risk value does not represent what actually happened but only stands for the relative relationships between districts. Figure 4 depicts the 5 × 5 semiquantitative matrix used for the risk assessment. Based on the matrix, the classifications of the risk assessment can also be graded, as shown in Table 1. Different risk values can be attributed to different degrees of risk.

Indices of Model Performance
To ascertain the simulated discharge values using the HEC-HMS model and the simulated water level using the FLO-2D model, three indices, flood peak error percentage (EYp), flood peak arrival time error (ETp), and Nash-Sutcliffe efficiency coefficient (NSE), were employed for model performance. These indices can be calculated as follows: (11) where N denotes the total number of observational/simulation data, The NSE value ranges from negative infinity to 1, and when the NSE is close to 1, it indicates that the model is of good quality and that the model has high reliability. If the NSE value is close to 0, the simulation result is close to the average level of the observed values; that is, the overall result is credible, but the error in the simulation process is large. If the NSE value is far less than 0, the model is not credible.

Results
To ascertain the availability of the models, two typhoon events, Typhoon Kongrey, which occurred from August 27 to September 2, 2013, and Typhoon Megi, which hit central Taiwan in the period of September 25-28, 2016, were adopted for the calibration and validation of the models, respectively. The validated models, HEC-HMS and FLO-

Indices of Model Performance
To ascertain the simulated discharge values using the HEC-HMS model and the simulated water level using the FLO-2D model, three indices, flood peak error percentage (EY p ), flood peak arrival time error (ET p ), and Nash-Sutcliffe efficiency coefficient (NSE), were employed for model performance. These indices can be calculated as follows: where N denotes the total number of observational/simulation data, (Y p ) obs and (Y p ) sim represent the observed and simulated peak discharge/water level, respectively, (T p ) obs and (T p ) sim indicate the observed and simulated peak arrival time for the discharge/water level, respectively, Y obs (t i ) and Y sim (t i ) denote the observed and simulated discharge/water level at time t i , respectively, and Y obs is the average value of the observations (= 1 The NSE value ranges from negative infinity to 1, and when the NSE is close to 1, it indicates that the model is of good quality and that the model has high reliability. If the NSE value is close to 0, the simulation result is close to the average level of the observed values; that is, the overall result is credible, but the error in the simulation process is large. If the NSE value is far less than 0, the model is not credible.

Results
To ascertain the availability of the models, two typhoon events, Typhoon Kongrey, which occurred from August 27 to 2 September 2013, and Typhoon Megi, which hit central Taiwan in the period of 25-28 September 2016, were adopted for the calibration and validation of the models, respectively. The validated models, HEC-HMS and FLO-2D, were then applied to predict inundation extents and depths, and these values were used to evaluate the hazard indices for different return periods.

HEC-HMS Model Calibration and Validation
To ensure the correctness of the rainfall runoff simulated by the HEC-HMS model, calibrating and validating the HEC-HMS model are important procedures. The observed discharges at the Wulilin Bridge flow station were compared with the simulated results. Figure 5a displays the comparison between the observed and simulated discharges at the Wulilin Bridge flow station during Typhoon Kongrey (2013); this comparison was used for the model calibration. The figure shows that the simulated results reproduced temporal variations in the discharge. Twin observed and simulated peak discharges were selected for the analysis of statistical errors. The EY p values for peak discharge 1 and peak discharge 2 were 0.6% and 19.3%, respectively, while the ET p values were 1 h and 0 h, respectively. The NSE value was 0.89.
2D, were then applied to predict inundation extents and depths, and these values were used to evaluate the hazard indices for different return periods.

HEC-HMS Model Calibration and Validation
To ensure the correctness of the rainfall runoff simulated by the HEC-HMS model, calibrating and validating the HEC-HMS model are important procedures. The observed discharges at the Wulilin Bridge flow station were compared with the simulated results. Figure 5a displays the comparison between the observed and simulated discharges at the Wulilin Bridge flow station during Typhoon Kongrey (2013); this comparison was used for the model calibration. The figure shows that the simulated results reproduced temporal variations in the discharge. Twin observed and simulated peak discharges were selected for the analysis of statistical errors. The EYp values for peak discharge 1 and peak discharge 2 were 0.6% and 19.3%, respectively, while the ETp values were 1 hr and 0 hr, respectively. The NSE value was 0.89.
A comparison of the simulated and observed discharges during Typhoon Megi (2016), used for the model validation, is shown in Figure 5b. The figure indicates that the computed water levels faithfully reproduced the observations. Based on the analysis of statistical errors, the EYp and ETp values were 4.6% and 2 hrs, respectively. The NSE value was 0.91.

FLO-2D Model Calibration and Validation
The observed water levels at the Changrun Bridge flow station during Typhoon Kongrey (2013) and Typhoon Megi (2016) were used to calibrate and validate the FLO-2D model, respectively. A comparison of the simulated and observed water levels during Typhoon Kongrey, used for the model calibration, is presented in Figure 6a. The figure A comparison of the simulated and observed discharges during Typhoon Megi (2016), used for the model validation, is shown in Figure 5b. The figure indicates that the computed water levels faithfully reproduced the observations. Based on the analysis of statistical errors, the EY p and ET p values were 4.6% and 2 h, respectively. The NSE value was 0.91.

FLO-2D Model Calibration and Validation
The observed water levels at the Changrun Bridge flow station during Typhoon Kongrey (2013) and Typhoon Megi (2016) were used to calibrate and validate the FLO-2D model, respectively. A comparison of the simulated and observed water levels during Typhoon Kongrey, used for the model calibration, is presented in Figure 6a. The figure shows twin peaks in the water level as a result of runoff discharge. The statistical errors represented by the EY p values at peak water level 1 and peak water level 2 were 13.5% and 23.1%, respectively, while the ET p values were 1 h and 4 h, respectively. The simulated water level indicated an underestimation of the observed water level at peak water level 1, while the simulation overestimated the observation at peak water level 2. The NSE value was 0.69.
shows twin peaks in the water level as a result of runoff discharge. The statistical errors represented by the EYp values at peak water level 1 and peak water level 2 were 13.5% and 23.1%, respectively, while the ETp values were 1 hr and 4 hrs, respectively. The simulated water level indicated an underestimation of the observed water level at peak water level 1, while the simulation overestimated the observation at peak water level 2. The NSE value was 0.69. Figure 6b illustrates a comparison between the computed and observed water levels during Typhoon Megi, as used for the model validation. The computed water levels matched the tendencies of the observations. The EYp and ETp values were 4.7% and 0 hr, respectively, and the NSE value was 0.79. The model performance observed during the model validation was better than that observed during the model calibration.

Inundation Depths and Extents under Different Return Periods
To obtain rainfall information under different return periods, this study used historical rainfall data from eight rainfall stations in the catchment area from 2012 to 2017 for frequency analysis, applied the Horner equation to establish the intensity-durationfrequency curve [43], and finally employed the alternating group method to obtain the designated rainfall of each rainfall station under different return periods. The designated rainfall was then input into the validated HEC-HMS model to obtain the runoff discharges in the upstream catchment area; then, the runoff discharges served as inputs into the validated FLO-2D model, and inundation disasters under different return periods regarded as hazards were obtained. For example, the analysis result for the rainfall hyetograph under a 200-year return period is illustrated in Figure 7.  Figure 6b illustrates a comparison between the computed and observed water levels during Typhoon Megi, as used for the model validation. The computed water levels matched the tendencies of the observations. The EY p and ET p values were 4.7% and 0 h, respectively, and the NSE value was 0.79. The model performance observed during the model validation was better than that observed during the model calibration.

Inundation Depths and Extents under Different Return Periods
To obtain rainfall information under different return periods, this study used historical rainfall data from eight rainfall stations in the catchment area from 2012 to 2017 for frequency analysis, applied the Horner equation to establish the intensity-duration-frequency curve [43], and finally employed the alternating group method to obtain the designated rainfall of each rainfall station under different return periods. The designated rainfall was then input into the validated HEC-HMS model to obtain the runoff discharges in the upstream catchment area; then, the runoff discharges served as inputs into the validated FLO-2D model, and inundation disasters under different return periods regarded as hazards were obtained. For example, the analysis result for the rainfall hyetograph under a 200-year return period is illustrated in Figure 7. 6.37 m, respectively. The maximum flooding site is located in the Yanchao District, resulting from the low-lying area around the riverside in the middle reach. There are 9 villages in the 50-year rainfall return period, 11 villages in the 100-year rainfall return period, and 12 villages in the 200-year rainfall return period with inundation depths greater than 3 m. Because most flooding zones are located downstream of the areas of the river where the flooding depth exceeds 0.3 m, the possible impacts on hazards are not negligible.

Inundation Depths and Extents under Different Return Periods
The inundation depths and extents were converted to hazard maps according to 5 grades: below 0.3 m, 0.3 m to 0.5 m, 0.5 m to 1.0 m, 1.0 m to 3.0 m, and above 3.0 m (see Table 1). Figure 9 presents hazard maps under the 50-, 100-, and 200-year rainfall return periods. The figure shows that there are 7 villages in Qiaotou District, 2 villages in Gangshan District, and one village each in Dashe, Nanzi, and Yanchao Districts belonging to the high and extremely high hazard levels under different rainfall return periods.

Inundation Depths and Extents under Different Return Periods
The inundation depths and extents were converted to hazard maps according to 5 grades: below 0.3 m, 0.3 m to 0.5 m, 0.5 m to 1.0 m, 1.0 m to 3.0 m, and above 3.0 m (see Table 1). Figure 9 presents hazard maps under the 50-, 100-, and 200-year rainfall return periods. The figure shows that there are 7 villages in Qiaotou District, 2 villages in Gangshan District, and one village each in Dashe, Nanzi, and Yanchao Districts belonging to the high and extremely high hazard levels under different rainfall return periods. Under the 200-year rainfall return period, the hazard grades of 12 villages increase to the extremely high level.

Inundation Depths and Extents under Different Return Periods
The inundation depths and extents were converted to hazard maps according to 5 grades: below 0.3 m, 0.3 m to 0.5 m, 0.5 m to 1.0 m, 1.0 m to 3.0 m, and above 3.0 m (see Table 1). Figure 9 presents hazard maps under the 50-, 100-, and 200-year rainfall return periods. The figure shows that there are 7 villages in Qiaotou District, 2 villages in Gangshan District, and one village each in Dashe, Nanzi, and Yanchao Districts belonging to the high and extremely high hazard levels under different rainfall return periods. Under the 200-year rainfall return period, the hazard grades of 12 villages increase to the extremely high level.

Analysis of Social Vulnerability and the Social Vulnerability Map
This study used two questionnaires to calculate social vulnerability. The first was the fuzzy Delphi questionnaire, which was adopted to screen the vulnerability factors. Through this questionnaire, experts were asked to identify the main factors affecting social vulnerability in urban areas. The second was the AHP questionnaire. The main purpose of this questionnaire was to set the weight value of each selected social vulnerability factor. Because some factors had a greater degree of influence than others, it was necessary to determine the weight value of each factor through the AHP method.

Analysis of Social Vulnerability and the Social Vulnerability Map
This study used two questionnaires to calculate social vulnerability. The first was the fuzzy Delphi questionnaire, which was adopted to screen the vulnerability factors. Through this questionnaire, experts were asked to identify the main factors affecting social vulnerability in urban areas. The second was the AHP questionnaire. The main purpose of this questionnaire was to set the weight value of each selected social vulnerability factor.
Because some factors had a greater degree of influence than others, it was necessary to determine the weight value of each factor through the AHP method.

Process with Fuzzy Delphi (FD)
Four assessment elements, which included twenty social vulnerability factors, were selected for the expert questionnaire, as shown in Table 2. Through the 41 collected expert questionnaires and the FD screen, the social vulnerability factors with consensus values higher than 6.09 were selected. Therefore, ten main social vulnerability factors, including individuals over 65 years old, individuals under 14 years old, elderly individuals living alone, individuals with disabilities, rubber boats, mobile pumps, firefighters, refuge shelters, schools (kindergarten, elementary, and high schools), and nursing homes for elderly individuals (see Table 2), were selected for further APH analysis.

Process with the AHP
After completing the fuzzy Delphi questionnaire and selecting the main vulnerability factors, an AHP questionnaire was conducted to determine the weight of each vulnerability factor. A total of 40 questionnaires were completed and returned by experts. Table 3 shows the weights of the vulnerability factors and the cascading weights of the hierarchy. The sum of the cascading weights of the hierarchy is 1.0. The table also indicates that experts believe that being a firefighter is the most important social vulnerability factor. From the AHP method, the weight of each factor can be determined, and the SVIs of urban areas can be evaluated.

Social Vulnerability Map
According to Equations (4)-(6), the social vulnerability index (SVI) can be obtained. To normalize the factors according to Equation (7), the SVI value was converted to the I SVI value. Since five social vulnerability levels were graded based on their I SVI values (see Table 1), the social vulnerability map was built, as illustrated in Figure 10. The map shows that the social vulnerability in the Nanzi District is the highest, at level 5 (i.e., extremely high vulnerability). The main reason for this result is that there are more nursing homes for elderly individuals, individuals over 65 years old, and individuals under 14 years old in this region than in other regions, resulting in the highest social vulnerability. The social vulnerability values of the two villages in Qiaotou District are both level 3. This is due to there being more nursing homes for elderly individuals in this district than in other districts.

Risk Map
After completing the hazard map under different rainfall return periods and the social vulnerability map, a semiquantitative risk assessment approach was employed to construct a risk map. The map of risks under different rainfall return periods is displayed in Figure 11. It shows that six villages in Qiaotou District and one village in Nanzi District belong to medium risk and extremely high risk levels, respectively, under the 50-year rainfall return period. There are seven villages and one village in Qiaotou District that are at medium risk and high risk levels, respectively, under the 200-year rainfall return period. It can be noted that the hazard risk under the 200-year rainfall return period is more serious than that under the 50-year rainfall return period, resulting in increased risks under the 200-year rainfall return period. However, the risk levels of most villages in the districts consist of extremely low and low levels. Social vulnerability values express relative, comparative relationships and do not represent actual vulnerability values. A village with a higher degree of vulnerability means that the village comprises a disadvantaged group with a denser population or fewer social resources. Therefore, when a disaster occurs, more serious damage and losses would occur in a village with a higher degree of vulnerability.

Risk Map
After completing the hazard map under different rainfall return periods and the social vulnerability map, a semiquantitative risk assessment approach was employed to construct a risk map. The map of risks under different rainfall return periods is displayed in Figure 11. It shows that six villages in Qiaotou District and one village in Nanzi District belong to medium risk and extremely high risk levels, respectively, under the 50-year rainfall return period. There are seven villages and one village in Qiaotou District that are at medium risk and high risk levels, respectively, under the 200-year rainfall return period. It can be noted that the hazard risk under the 200-year rainfall return period is more serious than that under the 50-year rainfall return period, resulting in increased risks under the 200-year rainfall return period. However, the risk levels of most villages in the districts consist of extremely low and low levels.

Risk Map
After completing the hazard map under different rainfall return periods and the social vulnerability map, a semiquantitative risk assessment approach was employed to construct a risk map. The map of risks under different rainfall return periods is displayed in Figure 11. It shows that six villages in Qiaotou District and one village in Nanzi District belong to medium risk and extremely high risk levels, respectively, under the 50-year rainfall return period. There are seven villages and one village in Qiaotou District that are at medium risk and high risk levels, respectively, under the 200-year rainfall return period. It can be noted that the hazard risk under the 200-year rainfall return period is more serious than that under the 50-year rainfall return period, resulting in increased risks under the 200-year rainfall return period. However, the risk levels of most villages in the districts consist of extremely low and low levels.

Discussion
Flood-related disasters in urban areas are subject to several factors, such as urban development, construction of transportation facilities, extreme weather and climate change, and anthropogenic changes that alter vulnerability, hazards, and risks [59]. To control the impacts of disasters caused by extreme rainfall events, an adaptation strategy must be adopted to reduce these disasters. To cope with the impacts of extreme rainfall events on urban areas, this study formulated adaptation strategies for villages with high

Discussion
Flood-related disasters in urban areas are subject to several factors, such as urban development, construction of transportation facilities, extreme weather and climate change, and anthropogenic changes that alter vulnerability, hazards, and risks [59]. To control the impacts of disasters caused by extreme rainfall events, an adaptation strategy must be adopted to reduce these disasters. To cope with the impacts of extreme rainfall events on urban areas, this study formulated adaptation strategies for villages with high risks levels. Generally, two measures, including nonstructural measures and structural measures, are recommended to mitigate disasters [25,60].
According to the flood risk assessment, the social vulnerability and hazards under different rainfall return periods reach extremely high levels, resulting in extremely high risk levels yielded for a village in Nanzi District. This village possesses a large population over 65 years old and under 14 years old as well as many nursing home facilities for elderly individuals. This study suggests that when planning escape routes with nonstructural measures, complex routes should be avoided, and the routes should be as simple and clear as possible to avoid dispersion in the escape process. In addition, the inundation depths are relatively high at this village, so it is recommended that a flood detention pond is added with structural measures to relieve the flooding load of the channel and reduce the flood disasters.
The social vulnerability of a village in Qiaotou District is calculated to be at level 3 (i.e., medium level), and the risk in this village for each return period is level 3 for the 50-year return period, level 3 for the 100-year return period, and level 4 for the 200-year return period. The number of mobile pumps in this area is lower than those of others area; therefore, this study proposes that mobile pumps are purchased in this area to reduce the damage caused by flooding.
The method used to reduce the social vulnerability index is a nonstructural measure, while the method employed to reduce hazards is a structural measure. To reduce the risk of flooding, it is necessary to first examine the sources of high risk from social vulnerability or hazards, as measures taken to reduce the risk must be adapted to the local conditions. In any case, the social vulnerability map, the hazard map, and the risk map can provide decision-makers with a reference that can be considered for risk management.
Most studies that have discussed hazard assessments applied a combination of hydrologic and hydrodynamic models for flood simulations [27,29,54,61]. Additionally, there are different research methods used to discuss vulnerability [23,33,62]. Further, two expert questionnaires were used to determine the main social vulnerability factors and the weights of the social vulnerability factors in the study area. This kind of research method is rarely proposed or discussed, but is specific and feasible, and the methodology can be applied to other urban areas affected by floods.
The limitation of this study is that because 20-m DTM data are available, the regional flooding simulation cannot display the flooded state of the street and because the data of the social vulnerability factor are not detailed enough, so the flood risk map can only be presented in the village.

Conclusions
Typhoon-induced extreme rainfall frequently results in flood damage that threatens the lives and property of people and causes economic losses in urban areas in Taiwan. Apart from structural measures that prevent flood damage, risk assessments play important roles in disaster management. To explore the flood risk map for disaster mitigation, the hypotheses of this study include the flood risk level in a village to be same, no engineering works executed to reduce hazard level, and social vulnerability for different rainfall return periods also to be same.
The flood risk assessment conducted in this study was established based on rainfallrunoff (HEC-HMS) and flooding simulation (FLO-2D) models applied for the hazard analysis and the fuzzy Delphi (FD) method and analytic hierarchy process (AHP) applied for the social vulnerability analysis. The HEC-HMS and FLO-2D models were calibrated and validated with observational data. The results indicated that the model performances during both the model calibration and validation were acceptable. The validated HEC-HMS and FLO-2D models were utilized to predict inundation depths and extents in urban areas under different rainfall return periods, and these factors were used to build hazard maps.
Furthermore, to explore social vulnerability, first, the fuzzy Delphi method was used to select the main social vulnerability factors, and then the AHP was adopted to determine the weight of each social vulnerability factor and to calculate the SVI to build a social vulnerability map.
By means of a semiquantitative risk matrix, the social vulnerability map and hazard map were combined into a risk map, and the risk map was divided into 5 levels to comprehend the risk status of urban areas in the Dianbao River watershed under different rainfall return periods.
Most studies did not clearly state that how the social vulnerability was selected for flood risk analysis. The main contribution of this study was that the hydrologicalhydrodynamic model, the fuzzy Delphi, and the analytic hierarchy process were proposed to build the flood risk map in urban areas of southern Taiwan. Based on the flood risk map, the nonstructural measures and structural measures were recommended to mitigate disasters according to the regional disaster characteristics.
Due to the impacts of climate change, urban flooding caused by extreme rainfall is worsening. In future research, we will explore a risk assessment of urban flood disasters caused by climate change.