Sea-level rise and flooding in coastal riverine flood plains

Abstract

This work presents a method for calculating the contributions of sea-level rise and urban growth to flood risk in coastal flood plains. The method consists of hydraulic/hydrological, urban growth and flood-damage quantification modules. The hydraulic/hydrological module estimates peak annual flows to generate flood stages impacted by sea-level rise within flood plains. A model for urban growth predicts patterns of urbanization within flood plains over the period 2010–2050. The flood-damage quantification module merges flood maps and urbanization predictions to calculate the expected annual flood damage (EAFD) for given scenarios of sea-level rise. The method is illustrated with an application to the Tijuana River of southern California, USA, and northwestern Mexico, where the EAFD is predicted to increase by over US$100 million because of sea-level rise of 0.25–1.0 m and urban growth by the year 2050. It is shown that urbanization plays a principal role in increasing the EAFD in the study area for the range of sea-level rise considered.

Editor Z.W. Kundzewicz

Citation Garcia, E.S. and Loáiciga, H.A., 2013. Sea-level rise and flooding in coastal riverine flood plains. Hydrological Sciences Journal, 59 (1), 204–220.

Résumé

Ce travail présente une méthode de calcul des contributions de l’élévation du niveau de la mer et de la croissance urbaine au risque d’inondation dans les plaines inondables du littoral. La méthode est composée de plusieurs modules : hydraulique / hydrologie, quantification de la croissance urbaine et quantification des dégâts dus aux inondations. Le module hydraulique / hydrologie estime les débits de pointe annuels pour générer les niveaux de crue touchés par l’élévation du niveau de la mer dans les plaines inondables. Un modèle de croissance urbaine prévoit les modes d’urbanisation dans les zones inondables sur la période 2010–2050. Le module de quantification des dommages dus aux inondations fusionne les cartes de zones inondables et les prévisions d’urbanisation pour calculer les dommages annuels attendus du fait des inondations (DAAI) pour des scénarios donnés d’élévation du niveau de la mer. La méthode est illustrée par une application à la rivière Tijuana, au Sud de la Californie, aux Etats-Unis, et au Nord-Ouest du Mexique, où les DAAI devraient augmenter de plus de 100 millions de dollars en raison d’une l’élévation du niveau de la mer allant de 0,25 à 1,0 m, et de la croissance urbaine d’ici à l’an 2050. On montre que l’urbanisation joue un rôle principal dans l’augmentation des DAAI dans la zone d’étude pour la gamme d’élévation du niveau de la mer considérée.

Key words:

• sea-level rise
• river hydraulics
• floods
• urban growth
• flood damage
• GIS

Mots clefs:

• élévation du niveau marin
• hydraulique fluviale
• inondations
• croissance urbaine
• dommages dus aux inondations
• SIG

1 INTRODUCTION

Flood-plain settlement is inherently risky. Yet, the occupation of available flat land within flood plains remains an enticing option for a growing population. Coastal areas, in particular, are densely populated and economically vibrant enclaves. Migration to coastal cities was an important urban trend through the 20th century (Small and Nicholls 2003, McGranahan et al. 2007) that continues unabated. Not surprisingly, urban disasters and environmental hot spots are frequent in coastal areas (Dilley et al. 2005). Population growth and the ensuing modification of land use are changes to be considered in coastal vulnerability assessments (Nicholls et al. 2008). Land-use change driven by population growth has the capability to alter streamflow patterns in a manner that may heighten flood risk. The increase of peak flows following the conversion of vegetated land to less pervious surfaces is a case in point (Loáiciga 2001, 2003a).

A second factor to be considered in the analysis of the flood risk in coastal riverine flood plains is sea-level rise caused by climate change. A rising sea level produces a higher downstream boundary condition on streamflow in a river’s estuary. A raised downstream river stage causes the flood stage to become higher than it would otherwise be upstream of the estuary and along the flood plain, thus raising the flood risk to life and property. Modern global sea-level rise is primarily attributed to thermal expansion and melting ice. The Intergovernmental Panel on Climate Change (IPCC) reported that the global mean sea level has been rising and is likely to continue to do so (IPCC 2007). The global average rate of increase was 1.8 ± 0.5 mm/year for 1961–2003 and 3.1 ± 0.7 mm/year for 1993–2003. Nicholls and Cazenave (2010) report a measured trend in mean sea-level rise equal to 3.3 mm/year in the period 1993–2009. The latter authors reviewed various predictions of sea-level rise by the end of the 21st century, with some being as little as 0.30 m or as high as 1.80 m over the year-2000 mean sea level. At San Diego, California, an area of interest in this work, sea level has been measured since 1906, and the mean sea-level trend is increasing at a rate of 2.1 ± 0.2 mm/year (sea-level monitoring station 9410170, NOAA 2009). Previous work has explored the potential impacts of sea-level rise on coastal resources at the global (e.g. McGranahan et al. 2007) and regional scales (e.g. Loáiciga et al. 2011).

This work presents a method to quantify the contributions of urban growth and sea-level rise to flood risk in coastal riverine flood plains at the regional scale. The flood risk is quantified in terms of the expected annual flood damage (EAFD), with all economic figures based on the year-2010 dollar value. The method is general in its applicability. It is illustrated with data for flood-prone areas along the lower Tijuana River in San Diego County, southern California. The method relies on modules for river stage calculation, urban-growth prediction and flood-damage quantification. Because of limitations on the long-range predictive skill of detailed urban growth by available models, the case study of this article is over a horizon ending in 2050. This article’s method of analysis is applied with three scenarios of total sea-level rise achieved with linear temporal growth by 2050: 0.25-, 0.5- and 1.0-m increments over the 2010 baseline sea level. The 0.25–0.5 m range is consistent with the median of sea-level rise predictions made for the 21st century (Nicholls and Cazenave 2010, Shepherd et al. 2012). The 1.0 m magnitude of sea-level rise corresponds to the upper range of sea-level rise predicted to take place during the 21st century (Nicholls and Cazenave 2010) and seems justified given recent estimates of accelerated polar ice melting (Shepherd et al. 2012). The results of this study can be used to identify areas particularly vulnerable to floods and to prioritize adaptive or remedial measures to cope with the heightened flood risk. Its novelty resides in coupling sea-level rise scenarios with hydrologic/hydraulic modelling and urban growth predictions in a bi-national coastal riverine flood plain to assess the change in the EAFD in the coming decades.

2 APPROACH TO FLOOD RISK ASSESSMENT: CONTRIBUTIONS OF SEA-LEVEL RISE AND URBAN GROWTH

The EAFD in coastal riverine flood plains is expressed as a function (F) of two variables: sea-level rise (SLR) and urbanization (U) or EAFD = F(SLR, U). Other variables, such as flood-control measures implemented in a river basin, could be important. Sea-level rise and urbanization are the key factors pertinent in this work. Let SLR₀ denote the baseline sea level from which future rise in sea level is measured, and U₀ is the baseline pattern and degree of urbanization. Baseline sea level and urbanization correspond to a common starting year, say, 2010. The baseline EAFD is denoted by F(SLR₀, U₀). The change in the EAFD (ΔF_SLR+U) caused by a change in sea-level rise (SLR = SLR₀ + ΔSLR) and a change in urbanization (U = U₀ + ΔU) is given by:

(1)

A change in EAFD caused by sea-level rise while urbanization remains at its baseline level, U₀, is given by:

(2)

Equation (2) allows quantification of the effect that sea-level rise alone has on EAFD, assuming that urbanization remains unchanged. It may also be of interest to quantify the effect that the growth of urbanization would have on EAFD if sea-level rise was to remain at its baseline level, SLR₀. The latter effect is expressed mathematically as:

(3)

The sum of the individual effects of sea-level rise (equation (2)) and the growth of urbanization (equation (3)) on the change of EAFD relative to its baseline level is equal to the combined effect of urbanization growth and sea-level rise on EAFD (embodied by equation (1)) only when the EAFD is linear on SLR and U, that is, when EAFD = F(SLR, U) = aSLR + bU + c, where a, b and c are constants. Otherwise, the sum of equations (2) and (3) can be either smaller or larger than the quantity represented by equation (1).

The function F is unknown in mathematical form. Therefore, one must approximate the formulas shown in equations (1)–(3) numerically, in this case without attempting to fit a linear equation to the EAFD involving the coefficients a, b and c stated above. To accomplish the numerical approximation, it is helpful to define a series of numerical scenarios, each associated with an EAFD caused by a specified combination of sea-level rise and urbanization. Thus, the baseline scenario is the EAFD₀ caused by baseline sea level SLR₀ and baseline urbanization U₀ or:

(4)

The calculation of EAFD₀ requires the implementation of a specialized hydrological technique to establish the severity of flooding associated with SLR₀. Next, one calculates the damage that the flooding associated with SLR₀ causes in a flood plain when there is a baseline degree of urbanization U₀. The flooding and damage variables are combined to render the respective baseline EAFD, EAFD₀. At the other extreme, one defines a scenario wherein the considered sea level and urbanization are realized at the end of the planning horizon, say by 2050. Denoting the future sea level and urbanization at the end of the planning horizon by SLR and U, respectively, the corresponding EAFD is then:

(5)

The calculation of the EAFD expressed by equation (5) requires a suitable future sea level to be set and the degree of urbanization by the end of the planning horizon to be predicted. The latter task involves the implementation of an urban growth model. The severity of flooding at the end of the planning horizon is established using the same hydrological technique employed to assess the baseline EAFD₀, as shown in equation (4), but imposing changed streamflow conditions created by a higher sea level. The end-of-planning-horizon flooding variables and the economic flood-damage characteristics of the predicted urbanization are evaluated (see equation (6) below) to arrive at the future EAFD in equation (5). The merging of hydrological analysis for flooding assessment and of urban-growth modelling are carried out for specified scenarios of sea-level rise and urbanization with which to construct the various changes in EAFD expressed by equations (1–3). The underpinning of the hydraulic/hydrologic analysis is presented next, followed by a description of urban-growth projections and flood-damage quantification.

3 HYDRAULIC/HYDROLOGIC ANALYSIS

In the context of EAFD estimation, the hydraulic/hydrologic analysis comprises two tasks: (a) determination of the flood-frequency function for the study area, that is the relationship between the magnitudes of (annual) peak flows and the probability of those flows and (b) determination of the relationship among peak flows and the flood stages (or river water–surface elevations) within the flood plain. The flood-frequency function may vary along a river reach, as discussed below. The function relating flood stage and peak flow is called a rating curve and is specific to a location on a stream within a flood plain. That is, the rating curve in general varies from location to location along a river reach even if the streamflow within the reach is constant. Knowledge of items (a) and (b) allows the derivation of the flood-stage vs probability function. In Section 5, we explain how to convert the flood-stage vs probability function into annualized flood damage once the urbanization characteristics have been determined.

The flood-frequency function expresses the relationship between the annual peak flow (Q) in a flood-plain stream and the (cumulative) probability (p) that Q be equal to or less than a specific value q of annual peak flow or p = P[Q ≤ q]. The return period (T, in years) corresponding to the flow q of probability p is equal to T = 1/(1 – p). In general, the flood-frequency function varies along a river reach. The annual peak flow of a specified return interval (or probability) typically decreases under natural conditions with increasing distance measured in the upstream direction from a river basin’s outlet. This is because the drainage area that contributes to the runoff decreases with increasing distance measured from the downstream boundary of the river reach. For this reason, when calculating the EAFD in a flood plain of a river reach, one chooses a reference river station at which to calculate the flood peak vs probability (flood-frequency) function. Taking, for the sake of argument, the 100-year annual peak flow (exceedence probability 1 – p = 0.01) at the reference river station, it follows that the latter peak flow defines a concomitant set of streamflow values throughout the river reach under study. This set of streamflow values, in turn, produces a profile of 100-year flood stage in the river reach and its flood plain, which has 100-year flood damage associated with it. It is in this sense that the flood-frequency functions, stage vs probability functions and flood damage vs probability functions must be interpreted for a river reach and its flood plain. In gauged basins (i.e. those with streamflow gauges), the flood-frequency function is derived from time series of annual peak flows by fitting a suitable probability density to the flow data (Interagency Advisory Committee on Water Data 1982, US National Research Council 1988). In ungauged basins, however, one common approach to the flood-frequency function relies on regional regression equations relating peak flows of specified return period to regionalized variables (such as mean annual precipitation, drainage area, topographic relief; see e.g. Benson 1962, 1964, Cruff and Rantz 1965,  Waananen and Crippen 1977). The regression approach is exemplified in the case study of this article.

The rating curve is developed by means of a river hydraulics model, which, as done in this work, employs a given downstream boundary condition for the flood stage and an associated value of peak flow to carry out a calculation of flood stage in the upstream direction of a river until reaching a user-specified upstream station. The public-domain Hydrologic Engineering Center River Analysis System (HEC-RAS) is a widely used river hydraulics numerical model in the United States (Warner et al. 2009, Brunner 2010). It was applied in this case study to calculate flood stage corresponding to floods of various return periods. The HEC-RAS is a one-dimensional (1D) flow routing model, which assumes that the left overbank and right overbank feature the same water level as that in the main channel. This is justified in relatively narrow flood plains with predominant 1D river flow (see case study below). A comparison of 1D and 2D river hydraulics models was presented by Horitt and Bates (2002).

The last step in the hydraulic/hydrologic analysis is to construct flood maps for the study area. This is accomplished by merging topographic information and the calculated water-surface elevations to create a flood-stage profile that is interpolated to create flood-plain extents or flood maps. The depth of inundation is determined within a flood zone by subtracting the DEM-generated terrain elevations from the calculated water-surface elevations. Interpolation of flood stage creates a flood-extent coverage for the river reach under analysis. Many flood-extent coverages are calculated as dictated by the number of peak-flow events (e.g. 2-, 5-, 10-, …, 100-year or larger peak flows) and the number of downstream flood-stage scenarios (e.g. corresponding to baseline sea level or to baseline plus an additional amount of sea-level rise).

Figure 1 shows a generic description of the procedure leading to the derivation of the flood stage (H) vs probability (p) function (Fig. 1(c)) based on the flood-frequency (peak flow vs probability, p versus q, Fig. 1(a)) and rating curve (flood stage vs peak flow, H versus q, Fig. 1(b)). First choose a pair (p, q) in Fig. 1(a). The peak flow q is then related to a flood stage H in Fig. 1(b). Finally, the pair (p, H) is located in Fig. 1(c). Repeating this procedure for various pairs (p, q) and (q, H) leads to the derivation of the flood stage (H) vs probability function (p). The latter function corresponds to a specified location in the stream, in the same manner as the flood-frequency function is location-specific on a river reach. A suitable reference river station is chosen to represent the flood stage vs probability function with the purpose of simplifying the graphical presentation of the H versus p relationship within a flood plain. In actuality, a river-hydraulics model yields a flood-stage profile from the downstream to the upstream boundary locations in a river reach and its flood plain for a chosen probability of annual peak. Thus, we must think of flood plain-wide flood-stage profiles vs probability relationships when conceptualizing the phenomena of flooding in a probabilistic sense.

Fig. 1 Derivation of the river stage (H) vs probability (p) function (part (c)) based on (a) the flood-frequency function (q versus p) and (b) the rating curve (H versus q). The flood peak (q) and stage (H) correspond to a designated reference river station. After Loáiciga (2003a). See text for details.

The case study of this article presents examples of the graphs shown in Fig. 1.

4 URBAN-GROWTH MODELLING

Urban-growth modelling generates predictions of urbanization in a populated coastal area vulnerable to river floods and sea-level rise over a period of analysis (until 2050). The modelling of urban growth starts with the baseline urbanization U₀ and ends with the end-of-horizon urbanization U. The urbanization scenarios (U₀ and U) are coupled with flooding scenarios induced by sea-level rise to calculate the EAFD associated with various urbanization/flooding combinations. The modelling of urban growth has a long history and is customarily linked with geographic information systems (GIS) technology to model complex populated areas (see work by Tobler (1970), and more recently, by Capello and Nijkamp (2009), for a review).

In this study, we adopted the SLEUTH model, formerly the Clarke Cellular Automaton Urban Growth Model (Clarke et al. 1997, Silva and Clarke 2002), to make projections of urban growth. The acronym SLEUTH comes from its six data inputs: (terrain) slope, land use, (land) exclusion, urban extent, transportation (in the form of a road network) and hillshade (used only for a background in the display of model graphics). The SLEUTH model is a cellular automata model of urban growth and land-use change. Cellular automata models have been used extensively in representing geographic processes, including urban and regional modelling (Agarwal et al. 2000). Regardless of the many differences between urban areas, it has been suggested that fundamental elements are common to each of them and that variation of these common elements results in the patterns of urban growth that are particular to any city (Silva 2004). The SLEUTH model captures the dynamics of urban growth and land-use change by undergoing calibration based on topographic data, historic trends of urbanization and road networks from previous time periods. Topographic information is expressed by terrain slope. It is known that increasing terrain slope is synonymous with decreasing urbanization. The SLEUTH model allows the specification of exclusion layers (or maps), which impose spatial constraints on urban growth, such as disallowing urban spread into water bodies.

Previous periods of urban extent and transportation-network information are used to calibrate the SLEUTH model, thus preparing it for making projections of a region’s future pattern of urban growth. The SLEUTH model’s calibration process conducts Monte Carlo trials with a set of model coefficients that are improved by comparing sequential simulations of urbanization in an area with historic urbanization in that area, until an acceptable match between the simulated and historical urbanization rates is achieved. There are five coefficients that affect the influence of the SLEUTH model’s growth rules, which determine the amount of spontaneous growth, new spreading centres, the rate of edge growth in urban areas, and to what extent road presence influences growth. The primary functions of the coefficient values are as follows:

1. Dispersion coefficient: weights the likelihood of spontaneous or road-influenced growth

2. Breed coefficient: weights the likelihood that a “spontaneous growth” becomes a “spreading” centre too

3. Spread coefficient: determines the likelihood that a pixel within a spreading centre generates additional urban pixels

4. Slope coefficient: accounts for a threshold slope at which new construction can no longer occur

5. Road gravity coefficient: the weight that proximity of a road lends to growth, based on distance.

Upon calibration, the SLEUTH model produces a set of coefficients, which take integer values in the range [1,100] that best simulate historical growth for a region. This set of coefficients, in conjunction with the initial conditions of land use in an area, generates projections of future urban growth. Such projections rely on the execution of urban-growth rules constrained by imposed exclusion scenarios. See Clarke and Gaydos (1998) for additional details on application of the SLEUTH model. The case study in this article illustrates the intricacies of implementing urban-growth projections with the SLEUTH model.

5 FLOOD-DAMAGE QUANTIFICATION

The calculation of EAFD in a vulnerable flood plain proceeds in two steps: (a) the flood damage, D, vs the stage, H, function is combined with the H vs probability, p, function to produce the D versus p function, denoted by D(p); and (b) the expected EAFD within the flood plain is given by:

(6)

Figure 2 illustrates in generic form the calculation of the EAFD: choose a pair (D, H) (Fig. 2(a)); relate the stage H to its probability p (Fig. 2(b)); and the pair (D, p) is noted in Fig. 2(c). Repeat this procedure for several pairs (D, H) and (H, p) to develop the D versus p function. The EAFD given by the integral in equation (6) is equal to the area ABCEFA in Fig. 2(c). The flood damage in the (D, H) function shown in Fig. 2(a) is a flood plain-wide amount. The flood stage refers to a chosen reference location or river station within the flood plain that is part of the profile of flood plain-wide stage for a given flood probability. Flood damage is calculated in this study using the HAZUS Flood Loss Estimation approach developed by the Federal Emergency Management Agency (FEMA 2009). The river hydraulics analysis model, HEC-RAS, calculates the flood stage and water velocity along a specified reach of a flood plain. The year-2010 flood stage rises, by assumption, by an amount equal to the rise in mean sea level at the river–ocean boundary associated with a given scenario of sea-level rise. The rise in downstream (boundary) flood stage propagates in the upstream according to the laws of gradually varied, steady-state flow. Upon calculation of the flood-stage profile along the flood plain, water elevations are mapped on a topographic map of the flooded area to determine the depth of flooding everywhere within the flood plain. Maps of the depth of flooding and of water velocity are combined with land-use maps to calculate the flood damage using FEMA’s HAZUS Flood Loss Estimation approach. The implementation steps leading to the calculation of EAFD are found in the case study below.

Fig. 2 (a) The flood damage (D) vs flood stage (H) function is combined with (b) the stage (H) vs probability (p) function to produce (c) the flood damage vs probability function. The EAFD is equal to the area ABCEFA in (c). The damage (D) represents the total economic damages over the area of flooding under study. The flood stage (H) corresponds to a designated reference river station and is part of the flood plain-wide stage profile corresponding to a flood peak probability, after Loáiciga (2003a).

The change in EAFD from a baseline condition of urbanization, U₀, and flood potential associated with baseline sea level, SLR₀, to a condition of urbanization, U, and flood potential exerted by higher sea level, SLR, is expressed by:

(7)

Figure 3 represents the situation embodied by equation (7), where the baseline condition and the future condition are represented by the damage (D) vs probability (p) functions ABC and EFG, respectively. The change in the EAFD in this case is the area enclosed by the polygon ABCGFEA.

In the next section, a case study is presented with an application of the hydraulic/hydrologic analysis, projection of urbanization and EAFD calculation in a coastal riverine flood plain threatened by sea-level rise and recurrent high flows.

6 CASE STUDY: FLOOD-HAZARD ASSESSMENT IN A COASTAL RIVERINE FLOOPLAIN

6.1 Study region

The Tijuana River drains a 4451.6 km² (1718.8 mi²) watershed that straddles the California–Baja California section of the USA–Mexico border. Approximately two-thirds of the Tijuana River basin lies in Mexico with the remainder located in the United States, as shown in Fig. 4.

Fig. 3 The change in EAFD between a baseline condition (line ABC) and a future condition (line EFG) is the area enclosed by the polygon ABCGFEA. D denotes flood damage, and p is the probability of annual peak flow.

The Tijuana River basin’s annual precipitation is highly variable and averages 20.4 cm. The City of San Diego is located in the lower reach of the Tijuana River, bordering the Pacific Ocean. With a population of 1.4 million people (January 2010 census), it is one of North America’s largest and fastest-growing metropolitan regions. Approximately one million new residents are expected in the region by 2030. The majority of urban development has been concentrated along the coast (URS Corporation 2007). The highest point of the watershed is located at Mount Laguna, California (elevation 1900 m), and its lowest point is at sea level at the Tijuana River estuary, in San Diego County, California.

Figure 5(a) shows the subdivision of the Tijuana River basin into its 12 constituent sub-basins (the naming convention relies on the major rivers they contain), which together with their drainage areas are (1) Pine Valley (254.1 km²), (2) Upper Cottonwood (353.3 km²), (3) Lower Cottonwood (includes Río Alamar) (350.4 km²), (4) Campo Creek (430.2 km²), (5) Río Seco (516.1 km²), (6) Río Tijuana (246.0 km²), (7) El Florido (282.6 km²), (8) La Cienega (513.0 km²), (9) Las Palmas (510.8 km²), (10) Las Canoas (356.1 km²), (11) Las Calabazas (397.9 km²) and (12) El Beltrán (241.1 km²). Precipitation data needed in the hydrological analysis of peak flows were available for seven regions of the Tijuana River basin. Therefore, the original 12 sub-basins were regrouped into seven sub-basins matching the seven regions for which precipitation data were available. Specifically, the Pine Valley and Upper Cottonwood sub-basins were merged into the North of Barrett Dam sub-basin (607.4 km²); the El Florido and Las Palmas sub-basins were combined into the Florido-Palmas sub-basin (793.4 km²); and La Ciénaga, Las Calabazas, Las Canoas and El Beltrán were merged into the Upstream (U/S) Valle de Las Palmas sub-basin (1508.1 km²). The Campo Creek, Lower Cottonwood, Río Seco and Río Tijuana sub-basins remained unaltered. The seven sub-basin groupings for discharge calculations are shown in Fig. 5(b). Table 1 shows data for the seven grouped sub-basins.

Fig. 4 The Tijuana River basin. The upper-left panel shows the Lower Tijuana River basin. The main panel depicts the entire watershed encompassing parts of the USA and Mexico and shows the perimeters of the cities of San Diego (USA) and Tijuana (Mexico).

Table 1 Drainage areas and average annual rainfall in the grouped sub-basins of the Tijuana River basin.

Sub-basin name	Sub-basin drainage area (km^2)	Sub-basin average rainfall (cm/year)
Río Tijuana	246.0	21.5
Campo Creek	430.2	23.6
Río Seco	516.1	16.8
North of Barrett Dam	607.4	36.6
Lower Cottonwood	350.4	30.5
Florido-Palmas	793.4	16.5
Upstream Valle de las Palmas	1508.1	13.5
Total Tijuana River basin	4451.6	20.4

Previous flood control studies in the Tijuana River basin have been reported by Ocean Science and Engineering Inc. (1971), the US Army Engineer District (1976), West Consultants Inc. (1994), the International Boundary and Water Commission (1994), Winckell and Le Page (2003) and others.

6.2 Hydraulic and urbanization scenarios

The baseline downstream flood stage in the Tijuana River is controlled by tidal level at the river–ocean boundary. The baseline downstream flood stage was determined from records kept within the Tijuana Estuary (PERL 2001). The three adopted year-2050 future downstream flood-stage scenarios correspond to the following:

1. 0.25 m of sea-level rise plus the 1.7 m maximum recorded flood stage to reach 1.95 m of flood stage at the downstream boundary

2. 0.5 m of sea-level rise plus the 1.7 m maximum recorded flood stage to achieve a downstream flood stage equal to 2.2 m

3. 1.0 m of sea-level rise plus the 1.7 m maximum recorded flood stage for a total downstream flood stage equal to 2.7 m.

The baseline (2010) downstream flood stage corresponds to 1.7 m. This is the maximum water surface elevation in the Tijuana River recorded at the downstream boundary. Flood-risk assessments were carried out for all the combinations of downstream flood stage and urbanization scenarios (baseline, or year-2010, and future, or year-2050). The combination of high stage of the downstream boundary with the occurrence of peak annual floods results in calculated flood stages through the study reaches of the Lower Tijuana River that are higher than those that would occur for less severe conditions of downstream boundary and peak flows. However, one must realize that high tides occur twice every 24 hours and that high stages are likely to prevail for about one third of each day. Under these circumstances, it is highly likely that a flood peak will take place at a time of high downstream stage.

6.3 Calculated peak flows at the downstream boundary of the Tijuana River

Henderson (2000) compiled a data set consisting of 18 years of continuous monthly precipitation for 44 stations in and around the Tijuana River basin. The mean annual precipitation (in cm) in each of the seven sub-basins and the entire basin are given in Table 1. These precipitation values are required in the estimation of peak flows (q) of various return periods (2-, 5-, 10-, 25-, 50- and 100-year). Equations for peak flows were reported by Waananen and Crippen (1977) for the San Diego region and other regions of California. The peak-flow equations were derived by means of regional regressions relating mean annual precipitation (P) and basin area (A) to peak flows of specified return periods in gauged basins, to be applied in ungauged basins of similar hydroclimatic setting. Peak-flow regression equations for the San Diego region are:

(8)

(9)

(10)

(11)

(12)

(13)

where q_N (N = 2, 5, 10, 25, 50 and 100) denotes the n-year peak flow (in m³s^–1) in equations (8)–(13); A represents the drainage area in km²; and P is the mean annual precipitation (in cm) over the drainage area. The peak flows q₂, q₅, q₁₀, q₂₅, q₅₀ and q₁₀₀ were calculated with equations (8)–(13), wherein the total basin area and basin-wide mean annual precipitation were set to A = 4451.6 km² and P = 20.3 cm, respectively; the calculated peak flows (rounded to the nearest integer) were 27, 131, 272, 626, 1002 and 1465 m³ s^–1, respectively. These peak flows were input to the HEC-RAS river analysis model to generate flood-stage profiles with which to calculate the flood damage in the flood plain. The calculated flood-frequency function at the downstream boundary of the Tijuana River is shown in Fig. 6.

Fig. 5 (a) Sub-basins for original watershed delineation and (b) grouped sub-basins.

6.4 Calculated flood-stage profiles

Flood stages were calculated using the HEC-RAS model, which requires cross-sectional (geometric) information and specification of hydraulic characteristics and peak flows along the river reaches modelled. Figure 7 shows the placement of the cross-sections, as drawn in ArcMap for the lower Tijuana River basin.

Fig. 6 Flood-frequency function q versus p at the downstream boundary of the Tijuana River.

The cross-sections extend laterally from the main channel into the flood plain to encompass areas projected by the urban-growth model SLEUTH to have a probability of at least 80% of becoming urbanized by year 2050. Several of the areas likely to become urbanized have near circular shapes with open areas near the centre that reflect hilly conditions with steep slope that are not buildable. The full specification of the geometric, hydraulic, peak flow and downstream boundary flood stage allowed the computation of flood stages along the lower Tijuana River basin with HEC-RAS. The river cross-section at station 1132 (i.e. located 1132 m upstream from the river mouth) was chosen as a reference station to define flood stage, flood peak in the hydrologic/hydraulic/economic method illustrated in Fig. 2. Figure 7 highlights the placement of cross-section 1132, which is the seventh cross-section counting upstream from the lower boundary of the Tijuana River. The flood-frequency function shown in Fig. 6 is applicable at station 1132, given its proximity to the downstream boundary of the river reach under consideration.

Fig. 7 Cross-sections within the area of analysis in the HEC-RAS flood-stage calculation. The reference cross-section 1132 is highlighted. (It is the seventh cross-section counted upstream of the lower boundary of the Tijuana River.)

For each scenario of sea-level rise and specified peak flow event (say, the 100-year event or q₁₀₀), HEC-RAS calculates the flood stage at each of its cross sections, yielding a flood-stage profile along the reach of the Tijuana River under consideration. Flood stages calculated with HEC-RAS at reference cross-section 1132 for the three scenarios of downstream (boundary) flood stage are shown in Table 2.

Table 2 Water-surface elevations (flood stages, H) in m, at reference cross-section 1132 calculated with HEC-RAS for baseline sea level and three sea-level rise scenarios.

	Flow (m^3s^–1)	Sea-level rise scenario, H (m)
		Baseline	+0.25 m	+0.50 m	+1.0 m
q2	27	1.72	1.95	2.20	2.70
q5	131	1.96	2.10	2.24	2.71
q10	272	2.21	2.29	2.36	2.74
q25	626	2.62	2.64	2.68	2.89
q50	1002	2.93	2.95	2.96	3.08
q100	1465	3.23	3.24	3.25	3.31

6.5 Projections of urbanization in the study area

Input data for SLEUTH are raster grids of pixels that cover an area typically derived from historical maps. Generally, data for these pixels (i.e. urban or non-urban, road or non-road) are processed and assigned through a geographic information system (GIS). This article focuses on urban projections for San Diego County, California. Of these county-wide projections, those for the study area pertinent to this work, the Lower Tijuana River sub-basin, were considered for flood-risk assessment. The SLEUTH model’s parameters were calibrated and validated in a previous biodiversity study that projected urbanization in San Diego County (Syphard et al. 2011). The calibrated SLEUTH model from this previous study was implemented here.

Syphard et al. (2011) evaluated two data sets with the purpose of assessing their relative merits in explaining existing urban patterns in San Diego County when used as inputs to the SLEUTH model. These data sets were called SILVIS and SANDAG (Syphard et al. 2011). This article relies primarily on the former data set because results presented by Syphard et al. (2011) showed a higher rate and greater extent of urbanization near the coast, thus providing the “worst-case” future scenario for coastal settlement (Syphard et al. 2011). SILVIS is a national database in which housing density has been mapped in the continental USA with decadal frequency since 1940 (Hammer et al. 2007). Land use is input to SLEUTH as either urban or non-urban. A housing-density threshold of 28 residental or non-residential building units per km² was chosen to categorize land use as either urban or non-urban. Data to ascertain the status of land use as either urban or non-urban were obtained from the California Department of Forestry and Fire Protection (http://frap.cdf.ca.gov). Figure 8 shows year-2010 maps for San Diego County.

Fig. 8 Input initialization data for San Diego County in year 2010 for the SLEUTH model. Maps of (a) percent terrain slope; (b) urban land cover; and (c) road network.

The urban land cover in Fig. 8(b) maps areas that exceed the urban density threshold of 28 units/km² (the total county urban extent was 11.5% of the landscape). Historical data on San Diego County’s road network from 1960, 1993 and 1997 were used in the calibration of the SLEUTH model. The most recent data, from 2000, were used in initializing SLEUTH’s urban projection. The year-2000 data were obtained from the San Diego Association of Governments (SANDAG). Calibration of the SLEUTH model based on the SILVIS data for San Diego County was accomplished by Syphard et al. (2011), resulting in the following SLEUTH coefficients: diffusion coefficient: 1; breed coefficient: 1; spread: 95; slope resistance: 1; and road gravity: 6. Goodness-of-fit evaluation of the SLEUTH simulations to historical data created with these coefficients yielded an optimal SLEUTH metric (OSM, detailed by Dietzel and Clarke 2007) equal to 0.665. The realized OSM metric is a measure of the goodness-of-fit between model predictions and observed urbanization. The calibrated SLEUTH coefficients imply that urban growth over time is strongly controlled by the rate of expansion from existing urban settlements, that is, by the spread coefficient.

The SLEUTH model simulated urbanization in San Diego County, resulting in projections for 2050. Two land-exclusion scenarios, presented in Fig. 9, constrained the simulations of urbanization (after Syphard et al. 2011). The first, and least restrictive, Scenario 1 (Fig. 9(a)) allows urbanization to occur on any privately owned land that is not excluded from future development. Scenario 1 excluded 0.102 km² of public lands from future development. Scenario 2 (Fig. 9(b)) allows urbanization in all public and private lands and reserves that are not excluded from future development. Scenario 2 excluded 8.14 km² of public and private lands from future development.

Fig. 9 The two land-exclusion scenarios used with the SLEUTH model in the Lower Tijuana River sub-basin: (a) Scenario 1 allows urbanization to occur on any privately owned land that is not excluded from future development and (b) Scenario 2 allows urbanization in all public and private lands and reserves that are not excluded from future development. The land covered by scenarios 1 and 2 represents exclusion conditions that remain constant up until year 2050.

The land excluded by scenarios 1 and 2 remains unchanged up to 2050 in the SLEUTH model simulations. The exclusion scenarios were based on: (a) land management data from the California Department of Forestry and Fire Protection (http://frap.cdf.ca.gov); (b) SANDAG GIS information showing land available for future development; and (c) lands protected from future development based on the San Diego County’s Multiple Species Conservation Program (MSCP) (Greer 2004). The first and second scenarios excluded, respectively, 23% and 49% of land in San Diego County from future urbanization or 0.04% and 3.3% of the approx. 246 km² Lower Tijuana River basin. Existing urban areas were classified according to a modified Anderson Level I Land Use Classification system (McGinnis 2001). The classification included residential, commercial, industrial, transportation and mixed urban/built-up land categories (Anderson et al. 1976, Syphard et al. 2011). The existing urban area within the Lower Tijuana River basin totalled approx. 5.4 km².

The calibrated SLEUTH model executes its growth rules starting with baseline land-use conditions and is constrained by the exclusion scenarios to generate projections of urban growth. The SLEUTH model projections calculated the probability of urbanization of land by 2050, where designation of land as urbanized requires that it have at least 28 units/km² in the data cells that are not explicitly excluded from development. The SLEUTH-calculated and rasterized probabilities were imported into a GIS, and the rasterized cells were categorized into four classes:

1. existing urban areas;

2. zero probability of urbanization (zero chance of having at least 28 unit/km² by year 2050);

3. low probability of urbanization (0–80% chance of urbanization); and

4. high probability of urbanization (80–100% chance of urbanization).

Figure 10 shows the areas calculated to have a probability ≥80% of becoming urbanized by 2050 for the two exclusion scenarios 1 and 2, which in the Lower Tijuana River sub-basin were 9.8 and 5.6 km², respectively, or approx. 4% and 2% of the area, respectively. It is seen in Fig. 10 that exclusion scenario 1 results in a greater land area projected to become urbanized and more vulnerable to flood damage by 2050 within the Tijuana River basin. For this reason, scenario 1, the worst-case Scenario, was chosen as the urbanization projection utilized in this work’s calculation of the EAFD for various sea-level scenarios.

Fig. 10 Areas calculated by the SLEUTH model to have a probability equal to or higher than 80% of becoming urbanized by year 2050; (a) and (b) show the areas associated with exclusion scenarios 1 and 2, respectively.

6.6 Calculation of EAFD

Residential and non-residential properties were assigned US$ values as presented in the City of San Diego’s vulnerability analysis completed for flood mitigation work (URS Corporation 2007). The FEMA hazard identification software, HAZUS, was used to determine replacement values for the San Diego County infrastructure threatened by flooding (all figures at year-2010 US$ values). The dollar value varies based on property attributes, such as type of structure, private vs public ownership and whether residential vs non-residential. For example, residential and non-residential buildings were valued at US$240 573 and US$2 528 760, respectively. The property–parcel data set was obtained from cadastral maps available through the San Diego Geographic Information Source (SanGIS) (SanGIS Data Warehouse 2010). Economic damages were determined only for properties and do not attempt to estimate the cost of loss of life. Figure 11(a) and (b) shows the 2-year and 100-year flood zones, respectively, associated with the 0.5- and 1.0-m sea-level rise scenarios (flood zones for the 0.25-m sea-level rise case were not included to avoid excessive detail in Fig. 11). Figure 11 also depicts (in red) the areas with at least 80% probability of becoming urbanized by 2050. Comparison of the 2-year and 100-year flood zones demonstrates the larger spatial extent of flooding as the peak flow increases. Flood-zone maps depict geographically the flood stage (H) generated throughout the flood plain by a peak flow of probability p.

Fig. 11 Flood-zone maps for (a) 2-year peak flow with 0.5- and 1.0-m sea-level rise; and (b) 100-year peak flow with 0.5- and 1.0-m sea-level rise. The areas of land with at least 80% probability of becoming urbanized by 2050 are shown in (a) and (b) in red.

The flood-zone maps were overlain with the property-parcel layers using GIS to develop the D versus p functions for each of the sea-level scenarios (baseline sea level in 2010 and 0.25-, 0.5- and 1.0-m sea-level rise by 2050) and for year-2010 and year-2050 urbanization. Figure 12 shows (i) D versus p for 0.25-, 0.5- and 1.0-m sea-level rise in conjunction with year-2050 projected urbanization and (ii) D versus p for baseline sea level with year-2010 urbanization. Notice that flood damages for flow with a return period of 2 years (p = 0.50) or less are the same for all of sea-level rise scenarios.

Fig. 12 Flood damage (D, in year-2010 US$) vs probability (p) functions in the Lower Tijuana River calculated for baseline (Base level) sea level and 0.25-, 0.5- and 1.0-m sea-level rise. The baseline D versus p graph corresponds to year-2010 urbanization, while those for sea-level rise correspond to year-2050 urbanization.

The results are summarized in Table 3, which indicates that there would be an increase in EAFD of US$245, 305 or 354 million in the Lower Tijuana River basin between 2010 and 2050 if the sea level were to rise by 0.25, 0.5 or 1.0 m, respectively, and if urbanization as projected by the SLEUTH model materializes by 2050. These increases are determined with respect to the baseline EAFD in 2010. These are high-cost losses expected if urbanization occurs as envisioned in our analysis for the range of sea-level rise considered. In addition, it is seen in Table 3 that if urbanization were to remain as it was in 2010, the three sea-level rise scenarios would produce much smaller increases in EAFD (US$39, 52 and 61 million). The implication is that the EAFD increases much more due to urbanization than due to sea-level rise, for the given range of sea-level rise and predicted urbanization.

Table 3 Change in the expected annual flood damage (ΔEAFD) associated with three sea-level rise scenarios and two states of urbanization, relative to the 2010 EAFD (amounts in US$ million at 2010 value).

State of urbanization	Sea-level rise
	0.25 m	0.50 m	1.0 m
2050 urbanization	245	305	354
2010 urbanization	39	52	61

7 CONCLUSIONS

This work introduced a method to assess the economic cost of bi-national riverine flooding aggravated by sea-level rise and urbanization. Our results demonstrate higher costs associated with sea-level rise in a rapidly growing coastal area, the Lower Tijuana River basin. By incorporating urban growth in the assessment of future EAFDs, the study has yielded a comprehensive estimate of flood damages likely to be experienced by 2050. Specifically, our results have demonstrated that future scenarios of 0.25–1.00 m of sea-level rise could result in increases of EAFD in excess of US$100 million for year-2050 urbanization in the Lower Tijuana River basin. Urban growth model projections have been shown to be essential in the assessment of location-specific flood risk. This is a novel feature in this type of flood-assessment study in areas threatened by sea-level rise and undergoing rapid urbanization. Our results concerning the EAFD in the Lower Tijuana River basin indicate that it increases much more due to urbanization than it does due to sea-level rise, for the given range of sea-level rise and predicted urbanization. The preponderance of direct human impact in this case has a noteworthy resemblance to the findings of the effects of groundwater extraction and climate change on aquifer storage in the Edwards Aquifer of Texas (Loáiciga 2003b). The latter study reported that groundwater extraction was the primary factor causing aquifer overdraft when compared with changes in aquifer recharge caused by changing climate.

The choice of sea-level rise scenarios used in this study is consistent with the most recent reports of thermal oceanic expansion and polar ice melting. Sea-level rise of between 0.25 and 0.50 m is realistic and likely to occur, while the 1.0 m rise is on the high side of predictions made by sea-level rise experts. It is acknowledged that predictions of urban growth are always shrouded in uncertainty due to the complexity of social processes governing land-use changes. The SLEUTH urban-growth model predictions, for this reason, are given in terms of a specified probability that an urban growth pattern might take place. Furthermore, our analysis took a conservative stance in choosing a worst-case scenario of land-use development. This acknowledges the erratic nature of future administrative decisions concerning the dedication of land to competing uses. Furthermore, the hydraulic routing of peak flows assumed that a high level of the downstream flood stage exists simultaneously with the occurrence of peak flows of specified return intervals in the Tijuana River. This choice of assumptions used in the calculation of future EAFD produced damages that are higher than those that would result if less severe conditions had been imposed. This choice of assumptions represents the authors’ acknowledgement that there exist urbanization and hydraulic/hydrologic uncertainties in the modelling of future EAFDs, and that, therefore, it is wise to assume worst-case scenarios in estimating future flood impacts.

The economic analysis of flood damage applied in this work accounted for direct damage sustained by structures, providing a cost approximation that does not factor in future property appreciation and other, difficult-to-quantify, second-order phenomena. Current urban extents and their projections were examined in conjunction with parcel-property information and aerial photography to determine specific land-use type for dollar assignment of damages. The flood-damage dollar estimates assigned to the evaluated structures were based on county-wide averages, which is an approximation to site-specific valuations.

Overall, this article presented a method for identifying areas simultaneously vulnerable to sea-level rise and undergoing urbanization. GIS provided a platform with which to jointly examine urbanization projections and varying flood-plain extents. This work’s calculation of economic damages with the purpose of comparing possible future land-use scenarios holds definite potential applicability in finding adaptive strategies to cope with sea-level rise, urban growth and flood risk in coastal flood plains.