**1. Introduction**

Urbanisation is a growing concern in the present world. With over 55% of the world's population living in urban areas [1], urban expansion has dramatically influenced the change in urban land use. Land use change due to urbanisation results in more impervious surfaces that have considerable impacts on urban hydrology [2]. Urban expansion leads to large impervious surfaces that reduces rainwater infiltration, generating high surface runoff and peak flow, thus increasing the risks of urban flooding and waterlogging [3]. Urban development is one of the major causes of urban pluvial flooding, aggravated by poor urban drainage systems that become more severe with increasing frequency and

**Citation:** Shrestha, S.; Cui, S.; Xu, L.; Wang, L.; Manandhar, B.; Ding, S. Impact of Land Use Change Due to Urbanisation on Surface Runoff Using GIS-Based SCS–CN Method: A Case Study of Xiamen City, China. *Land* **2021**, *10*, 839. https:// doi.org/10.3390/land10080839

Academic Editors: Matej Vojtek, Andrea Petroselli and RaffaelePelorosso

Received: 14 July 2021 Accepted: 9 August 2021 Published: 11 August 2021

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

**Copyright:** © 2021 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

magnitude of extreme precipitations as a consequence of climate change, affecting the urban population, infrastructure, and economy [4–6].

The effect of land use change on surface runoff depends on many factors such as spatial characteristics of land use, rainfall, and soil [7]. The effect is most prominent in areas undergoing rapid urbanisation where natural vegetation has been removed and replaced by constructed land [8]. Removal of trees and vegetation and construction of commercial and residential buildings, streets, and parking lots alters the water balance of an area by changing the balance between rainfall, infiltration, evapotranspiration, and runoff [9,10]. Constructed areas such as concrete buildings, roofs, paved roads, and sidewalks lead to an increase in impervious surfaces, which reduces the time of concentration for runoff, leading to higher storm runoff and peak discharges [11,12]. As a result, flooding and waterlogging are common issues in highly urbanised areas [13].

The effect on hydrological processes, particularly surface runoff characteristics, of land use change due to urbanisation has been studied by many researchers [14–16]. Specifically, Hu, Fan, and Zhang [2] investigated the impact of land use change on the distribution of surface runoff in the highly urbanised city of Beijing. Results indicated that the change in surface runoff is strongly associated with change in impervious areas. Chen et al. [17] showed that non-uniform urban expansion and intensification are the major driving forces for changing surface runoff. The studies mentioned above indicate that changes in urban land use with more impervious surfaces reduce the ability to intercept rainfall and exacerbate surface runoff process.

The SCS–CN method, developed by the Natural Resources Conservation Service (NRCS), U.S. Department of Agriculture (USDA), is one of the most widely used methods in many hydrological studies, particularly for estimating surface runoff, accounting several factors such as soil, land use treatment, and topographical features, and incorporating these factors into a single CN parameter [18–21]. However, the conventional method is more time consuming and error prone. Remote sensing (RS) and a geographic information system (GIS) can be used effectively to manage spatial and non-spatial databases that represent the hydrologic characteristics of the watershed [22]. Therefore, the application of GIS and remote sensing with hydrological models yield results with high reliability and accuracy over conventional methods [23]. Many hydrological studies are being conducted with the integration of GIS and RS techniques for simulating surface runoff. For instance, Shadeed and Almasri [24] demonstrated that the integration of GIS with SCS–CN for simulating runoff volume proved to be an effective tool in arid to semi-arid catchments [24]. Similarly, Liu and Li [25] used a GIS-based SCS model for runoff generation in Loess plateau of China and found that the calculated and observed runoff processes were well correlated.

Recently, the GIS-based SCS–CN method has been used widely to study the impacts of land use change and urban growth on surface runoff [2,14,26]. For example, Jahan et al. [27] investigated the impacts of suburban growth on surface runoff using GIS-based SCS–CN, and concluded that the integrated approach is useful for land use change detection and analysis of impact of suburban growth on surface runoff. Furthermore, Vojtek and Vojtekova [14] applied the SCS–CN method in a GIS environment to study the impact of land use change on surface runoff. Results indicated that the method is suitable for spatial modelling of runoff characteristics of a basin where few data are available. Therefore, GIS and remote sensing in the SCS–CN model is a more efficient and cost-effective method that produces adequate results without using complex data in analysing the impact of land use change and urban growth [27–30]. In this study, the GIS and SCS–CN method is integrated to analyse the effect of land use change, and examine the spatial-temporal variation in surface runoff at an urban scale.

Firstly, the study assessed land use change in Xiamen city from 1980 to 2015 by analysing land use images in GIS environment. Then, surface runoff depth was estimated using SCS–CN method. We used GIS-based SCS–CN method to evaluate the effect of land use change on surface runoff. The main objective of this research is to identify the characteristics of land use change from 1980 to 2015 and analyse the impact of land use

change due to urbanisation on the temporal and spatial distribution of surface runoff under the rainfall return periods of 5, 10, 20, and 50 years.

#### **2. Materials and Methods**

*2.1. Study Area*

Xiamen City is located on the west coast of the Taiwan Strait, composed of the mainland area along Xiamen Bay, Xiamen Island, Gulangyu, and other islands, including 6 administrative districts: Huli, Siming, Haicang, Jimei, Tong'an, and Xiang'an, as shown in Figure 1. The city has a total land area of 1699 km<sup>2</sup> and a sea area of 324 km2. The study area is predominantly flat, low relief, with medium-low mountains, plains, and tidal flats. The slope of the terrain descends from northwest to southeast. The northwest part is mountainous, with the highest elevation of 1175 m above sea level, located on Yunding mountain [31]. Xiamen has a subtropical monsoon climate, with humid and mild climates throughout the year. According to the *Xiamen Statistical Yearbook*, annual average temperature and rainfall are approximately 21 ◦C and 1200 mm, respectively, with most rainfall in June, July, and August, accounting for more than half of the annual rainfall. Xiamen is one of the special economic zones in China and has experienced a rapid economic development and urbanisation [31]. According to Xiamen Municipal Bureau of Statistics [32], urban built-up area of the city has expanded from 38.5 km<sup>2</sup> in 1985 to 348.23 km<sup>2</sup> in 2017. The permanent population of the city was 4.29 million in 2019 [32]. As of 2019, the urbanisation rate of the population in Xiamen is 89.2% [33]. Xiamen has experienced several natural disasters, especially flooding and waterlogging induced by sea-level rise and storm surges, which is further intensified by rapid urban growth in past decades. Hence, there is an urgen<sup>t</sup> need for mitigation of the potential hazards induced by urban development and climate change.

**Figure 1.** Location of study area with administrative boundary and catchment divisions.

Figure 2 shows the overall methodology of the study conducted for analysing impact of land use change on surface runoff with the application of GIS and SCS–CN method. The GIS environment was used to build and intersect land use and HGS shapefiles. The land use and soil complex were used to obtain the weighted CN. The CN value of each polygon was estimated using the USDA table. Finally, raster calculator in the GIS was used to calculate runoff depth from CN values for the four rainfall return periods. The methodology was applied to four different years of 1980, 1990, 2005, and 2015, respectively.

**Figure 2.** Flowchart showing the methodology of the study.

#### *2.2. Data Source and Methods*

#### (a) Land use data

Multispectral satellite images from Landsat 3 MSS, Landsat 5 TM, and Landsat 7 ETM+ were obtained to create a land use map of the study area for 1980, 1990, 2005, and 2015. The images were downloaded from Centre for Earth Observation and Digital Earth (CEODE), Chinese Academy of Sciences, and the United States Geological Survey (USGS) [34]. Images were georeferenced to the UTM, Zone 50 North, WGS-84 projection, and Beijing 1954 coordinate systems. The spatial resolution is set to 30 m. The land use was classified according to the National Land Cover Data Sets (NLCD) of China generated by Liu et al. [35] in the construction of the *China 20th Century LUCC Spatio-temporal Platform*. The method of unsupervised classification with Iterative Self Organizing Data Analysis Technique Algorithm (ISODATA) was used for image classification. Based on the national land use data product by Geographical Information Monitoring Cloud Platform, the land use was initially classified as six first level categories and twenty-five second level categories, which were reclassified into eight classes for the study. High resolution satellite imagery from Google Earth was used for visual interpretation and accuracy assessment. The overall accuracy and kappa coefficient of classified land use types is more than 85%. The high overall accuracy and kappa coefficient suggests a good relationship between classified image and reference image. The detailed descriptions of land use classes are given in Table 1. ArcGIS 10.5 and ENVI 5.3 software were used to generate various layers and land use maps.


**Table 1.** Land use classes and their description.

#### (b) Soil data

Soil data were obtained from Harmonized World Soil Database (HWSD) and Food and Agriculture Organization (FAO). Soil information for Xiamen city was obtained from 1:1 million soil map of China with a resolution of 1 km (30 × 30 arc seconds) [36]. Hydrological soil groups are classified into A, B, C, and D on the basis of water transmission and infiltration rate of soil when the soil is thoroughly wetted [29]. Soil group A has a high infiltration rate and the lowest runoff potential. These soils consist of sands or gravels that have a high rate of water transmission. Group B and C soils have moderate runoff potential and infiltration rate with a slower rate of water transmission. Soil group D has the lowest infiltration rate and highest runoff potential composed of clay or clayey loam soils with a very slow rate of water transmission [29]. Based on soil texture and soil type, the study area is classified into three HSG types, B, C, and D, as can be seen in Figure 3, which mostly contribute to large surface runoff and less infiltration.

(c) Rainfall data

Long-time series rainfall data from 1985 to 2015, collected at the Xiamen meteorological stations, were acquired from the China Meteorological Data Service Centre of the China Meteorological Administration. The rainfall amount for the return periods of 5, 10, 20, and 50 years were obtained from maximum daily rainfall data for the hydrological analysis, as seen in Table 2. The rainfall for corresponding return periods was determined by using Log Pearson Type III distribution, which is the common probability distribution method used in China [37]. The study assumed that the climatic and soil conditions are constant. To investigate spatial heterogeneity, the study required uniformly distributed station data, which are not available. Therefore, for the temporal analysis of rainfall from the stations, the station with the most complete data was selected as a representative of the study area's rainfall. Furthermore, the study area is small, and changes in rainfall variation is considered insignificant.

(d) Digital Elevation Model (DEM)

A digital elevation model (DEM) with a 30 m resolution was downloaded from the USGS website. Using hydrology tools in ArcGIS 10.8, 19 catchments were obtained, as seen in Figure 1. The sinks in the DEM data were filled, water flow direction was estimated, and flow accumulation was set with a threshold of 10,000 and 60,000 to obtain the catchments. As the study area has short streams and not many large rivers, areas in the southwest do not have adequate streams to obtain a catchment. Hence, runoff obtained from these areas

was not significant, and only those areas with enough streams within the threshold are selected as catchments for runoff analysis.

**Figure 3.** Soil type and hydrological soil group in Xiamen.

**Table 2.** Rainfall depth for 24 h maximum daily rainfall in different return periods.


Source: Xiamen Meteorological Bureau.

#### *2.3. SCS–CN Method*

The SCS–CN method is the most commonly used empirical hydrological method developed by the NRCS, USDA, and is widely used to simulate runoff [19,29]. The curve number method has been successfully adopted in many ungauged watersheds and has expanded its scope of application in urbanised catchments and forested watersheds [29]. The surface runoff model uses the curve number approach of the US Soil Conservation Service [29], based on combinations of land use, hydrological soil group, and antecedent moisture condition (AMC) for the estimation of runoff. The amount of runoff was estimated using the SCS–CN method in presence of GIS and RS. The curve number is the most important factor in determining runoff via the SCS based method. The runoff of the soil and land use complex is represented by CN, which is a function of soil type, moisture conditions, and land use type [38].

The SCS–CN model is based on the water balance equation as shown in the Equations (1)–(3).

 −

$$P = Ia + F + Q \tag{1}$$

$$\frac{Q}{P - Ia} = \frac{F}{S} \tag{2}$$

$$Ia = \lambda \times \mathcal{S} \tag{3}$$

where *P* is the rainfall depth (mm), *Ia* is the initial abstraction of the rainfall (mm), *F* is infiltration, *Q* is surface runoff depth (mm), *S* is the potential maximum soil retention, and *λ* is abstraction coefficient that ranges between 0.0 to 0.2 or 0.05 for urbanised catchments [39]. The value of 0.2 as mentioned by Natural Resources Conservation Service (NRCS) was used in the study [29]. The runoff depth can be obtained for two conditions from the Equations (1) and (2):

For *P* > *Ia* and

$$Q = \frac{\left(P - Ia\right)^2}{P - Ia + S} \tag{4}$$

 If *P* < *Ia*, *Q* = 0 and *Q* from Equation (4) is expressed as follows:

$$Q = \frac{\left(P - 0.2S\right)^2}{P + 0.8S} \tag{5}$$

In Equation (5), *S* was obtained from the dimensionless parameter *CN*. *CN* is runoff curve number that ranges from 0 to 100.

$$S = \frac{25,400}{\text{CN}} - 254\tag{6}$$

In the SCS–CN method, curve number plays an important role in determining the surface runoff of an area, and its value depends on the corresponding soil type and AMC. AMC is antecedent moisture condition present in the soil at the beginning of the rainfall. Areas with a higher curve number represent higher runoff generated from the surface. In our study, we chose B, C, and D as three hydrological soil groups (HSG) present in the study area, and the soil moisture condition (AMC II) was set as moderate according to average runoff condition. A combined map of land use and HSG was generated by combining land use and soil maps in ArcGIS using overlay analysis. Then, the CN values were assigned for each polygon based on the information of land use and soil. The CN values for each land use type under AMC II is obtained from the TR-55 lookup Table 3 [29,40].

**Table 3.** CN numbers for corresponding land use types.


Area weighted CN was obtained to simulate the runoff of the whole area using the initial curve numbers from the table as in Equation (8). Combining the CN values of different land use and soil complex polygons, weighted CN was calculated for each catchment. The weighted CN is calculated by taking the sum of each CN value multiplied by its fraction of the total area of each land use type [41]. The Equation is given below:

$$CNw = \frac{\sum (CNi \times Ai)}{A} \tag{7}$$

where *CNw* is the weighted curve number; *CNi* is the curve number for each land use type; *Ai* is area of land use with respective curve number; and *A* is the total area of each land use type. Finally, surface runoff depth was estimated, and runoff coefficient, i.e., the ratio of runoff to rainfall, was calculated.

#### *2.4. Analysing Impact of Land Use Change on Surface Runoff*

Impact of land use change on surface runoff was analysed by comparing the difference in runoff variables. Runoff depth and runoff coefficient was used as two variables to assess land use change on surface runoff. Surface runoff from the catchments were obtained for the different land use conditions and runoff was calculated under different rainfall return periods. The difference of runoff and runoff coefficient was obtained for the land use period of 1980–1990, 1990–2005, and 2005–2015 by using Equations (8)–(10):

$$
\Delta Q = Qa - Qb \tag{8}
$$

$$
\Delta \ll = \frac{\Delta \bar{Q}}{P} \tag{9}
$$

$$
\Delta \beta = \frac{\Delta Q}{Qb} \times 100\% \tag{10}
$$

where *Qa* and *Qb* denotes surface runoff depths (mm) of the initial and final land use, respectively, *P* is the rainfall depth (mm), Δ*Q* is the change in runoff depth between two periods of land use conditions, Δ*α* is the absolute change in the runoff coefficient, and Δ*β* represent relative change in runoff. Land use change leads to an increase in surface runoff if the values of Δ*Q* and Δ*α* are positive in the above Equation (10).

A relationship between land use and surface runoff was determined by using Pearson's product–moment correlation coefficient. A positive and larger correlation coefficient suggests that the factor is more significant in the change in surface runoff.

#### *2.5. Validation of SCS–CN Model*

For analysing the performance of the model, the model was validated using observed flow and simulated runoff between 1981 and 2015. Four statistical indices, as shown in Table 4, were used for testing the goodness of fit.


**Table 4.** Statistical indices used for model validity.

## **3. Results**

*3.1. Land Use Change in the Study Area*

Land use types for four different years can be seen in Figure 4. The spatial distribution of land use types show that the area has experienced significant change, especially due to urban expansion. There was a remarkable expansion of built-up areas from 62.85 to 307.54 km<sup>2</sup> between 1980 and 2015. It was the largest gain, with a net increase of 244.69 km2. The largest net loss from 1980 to 2015 was observed in farmland and coastal wetlands with 218.2 km<sup>2</sup> and 46.65 km2, respectively. However, constructed land that includes built-up land and rural settlements underwent the largest net increase. The forest area reduced substantially with a change of 23.93 km2.

**Figure 4.** Land use maps of Xiamen in 1980, 1990, 2005, and 2015.

Similarly, grassland and farmland were replaced by built-up areas in the south and south-eastern parts of the mainland and major parts of the island. It can be observed that most coastal wetlands decreased from 1980 to 1990, and built-up areas occupied most areas. The land area of Xiamen City expanded outward by reclamation and construction along the coast after the national survey of 1985, which contributed to the expansion of land area and reduced wetlands in 1990 [42]. A considerable increase in water bodies indicate some farmlands being converted to reservoirs, aqua farms, and other constructed wetlands. Table A1 shows the conversion of the land use in the form of a change matrix for the period of 1980 to 2015. There was a major conversion from farmland to constructed land. The trend of land use change from 1980 to 2015 indicates that urban development has dominated the island city, resulting in vast areas of impervious surfaces.

In Table A2, the percent and area change of each land use type is depicted from 1980 to 2015. It can be observed that farmlands have been reduced by 14.01% and built-up area has increased by 15.7%. The transfer process of land use types in between the study years can be seen from Figure 5. It can be noticed that major portion of change in farmlands in 1990 and 2005 have been replaced by built-up land. Overall, constructed land increased from 9.12% in 1980 to 26.1% in 2015. This indicates that urban impervious areas have increased considerably in the last few decades of the study period.

**Figure 5.** Transfer process of land use type during 1980–2015.

#### *3.2. Spatial Distribution of Runoff in Different Years*

The spatial distribution of runoff depth for the highest rainfall return period of 50 years is shown in Figure 6. Surface runoff depth ranges from 176.28 mm to 329.16 mm. The area covered by built up land with higher CN depicts a higher amount of runoff in the period from 1980 to 2015. Xiamen Island, particularly dominated by constructed area, shows larger areas of high runoff. The urban areas with high CN are often the areas with high runoff value.

Under the land use conditions of 1980, 1990, 2005, and 2015, the average surface runoff depth and runoff coefficient show an increasing trend as shown in Figure 7. The average surface depth of the area under four different periods differ from 117.2 to 271.6 mm, and runoff coefficient fluctuated from 0.6 to 0.8. The calculated value of average surface runoff and runoff coefficient increases as the rainfall return period increase from 5 years to 50 years. It can be observed that surface runoff and runoff coefficient significantly increase from the year 1990 to 2015. The amount of runoff percent for different land use types can be noticed in Figure 8. A major portion of runoff is contributed by built-up land and rural settlements which are the major constructed land in the study area with the rise from 38.2% to 48.4%. Runoff contributed by built-up land alone in the study area increased from 14.2% to 27.9% from 1980 to 2015. Similarly, farmlands and coastal wetlands also contribute to the significant portion of surface runoff in the study area.

**Figure 6.** Surface runoff depth in Xiamen from 1980 to 2015.

**Figure 7.** Average surface runoff depth Q (mm) and runoff coefficient (RC) under the return period of 5, 10, 20, and 50 years under the land use conditions of 1980, 1990, 2005, and 2015.

**Figure 8.** Runoff depth across different land use types.

#### *3.3. Change in Surface Runoff in Different Land Use Conditions*

Increase in urbanised area gives rise to impervious surfaces, resulting in increased surface runoff. Table 5 depicts the change in surface runoff at different stages of land use conditions for the rainfall return periods of 5, 10, 20, and 50 years. Increase in amount of runoff (ΔQ) is highest during 1990–2005 and lowest during 1980–1990. The surface runoff coefficient change (Δα) during the periods is influenced by the respective change in average surface runoff under return periods of 5, 10, 20, and 50 years.

**Table 5.** Changes in surface runoff at different stages under rainfall return periods of 5, 10, 20, and 50 years.


#### *3.4. Relationship between Surface Runoff and Land Use*

Land use exhibits a high relationship with surface runoff with a statistically significant correlation (*p* > 0.05). As shown in Table 6, farmlands (−0.97), forestland (−0.96), and grassland (−0.97) contribute negatively to the surface runoff, which indicates that increase in these land use types contribute to a decrease in surface runoff. Increase in urban built-up land (0.98) and rural settlements (0.99) corresponded to an increase in average surface runoff, while the decrease in farmland, forestland, and grassland contribute to higher runoff. Furthermore, increase in coastal wetlands contribute to a decrease in runoff, and unused land contributed positively to average runoff. The results are consistent with the common knowledge that increase in urban constructed land causes an increase in surface runoff; however, increases in farmland, grassland, forestland, and wetlands lead to a reduction in surface runoff. Water bodies are considered to have less effect on runoff depth, hence the relationship is less significant.

**Table 6.** Relationship between land use and surface runoff under rainfall return period of 50 years. t—t-statistic, *p*—*p*-value, and R2—Pearson's correlation coefficient. Areas in square kilometres.


*3.5. Rainfall-Runoff Correlation Analysis*

The correlation analysis shown in Figure 9 indicates a strong linear relationship between the SCS–CN runoff and maximum daily rainfall, with a correlation coefficient of 0.99. The study results are comparable with the findings of Rawat and Singh [43] who found a good coefficient of determination (0.91) in a small study area using the SCS–CN model. The slope of the line determines the runoff coefficient, i.e., 0.84. The resulting findings are similar to Al. Ghobari et al. [44], who came to the conclusion that the SCS–CN model has a better simulation effect on study areas with a coefficient of runoff greater than 0.5 than those with a coefficient of runoff less than 0.5. This coefficient may provide valuable information on the extent of the basin response to runoff generation.

**Figure 9.** Relationship between rainfall and SCS–CN runoff.

#### *3.6. Validation of SCS–CN*

The SCS–CN model was validated using historical observations and simulated flow from 1980 to 2015, as shown in Figure 10. Some statistical efficiency criteria are used to perform evaluation of the validation results between simulated output and observed data which are percent bias (PBIAS), correlation coefficient (r), Nash–Sutcliffe efficiency (NSE) and volumetric efficiency (VE). These statistical indices indicate the goodness of fit between

simulated and observed data. The model successfully predicts the annual flow with the high accuracy as depicted by the indices. The PBIAS, r, NSE, and VE were −5.7, 0.82, 0.64, and 0.86, respectively. Although there is an underestimation of streamflow due to static land use, annual flow statistics indicate that there is a good relationship between observed and simulated streamflow. Hence, the model performance was satisfactory and responded well in simulation of runoff.

**Figure 10.** Observed and simulated annual flow from 1981 to 2015.
