Assessing Sponge Cities Performance at City Scale Using Remotely Sensed LULC Changes: Case Study Nanjing

: As a result of high-density urbanization and climate change, both the frequency and intensity of extreme urban rainfall are increasing. Drainage systems are not designed to cope with this increase, and as a result, ﬂoods are becoming more common in cities, particularly in the rapidly growing cities of China. To better cope with more frequent and severe urban ﬂooding and to improve the water quality of stormwater runoff, the Chinese government launched the national Sponge City Construction (SCC) program in 2014. The current SCC design standards and guidelines are based on static values (e.g., return periods, rainfall intensities, and volume capture ratio (VCR)). They do not fully acknowledge the large differences in climate conditions across the country and assume that the hydraulic conditions will not change over time. This stationary approach stems from the traditional engineering approach designed for grey infrastructure (following a “one size ﬁts all” approach). The purpose of this study was to develop a methodology to assess the VCR baseline (before construction in the pre-development stage) and changes in VCR (difference between the VCR of the pre- and post-development stage). The VCR of the post-development stage is one of the required indicators of the Assessment Standard for Sponge Cities Effects to evaluate SCC projects. In this study, the VCR was derived from remote-sensing-based land use land cover (LULC) change analysis, applying an unsupervised classiﬁcation algorithm on different Landsat images from 1985 to 2015. A visualization method (based upon Sankey chart, which depicts the ﬂows and their proportions of components) and a novel and practical partitioning method for built-up regions were developed to visualize and quantify the states and change ﬂows of LULC. On the basis of these ﬁndings, we proposed a new indicator, referred to as VCRa − L , in order to assess the changes in urban hydrology after SCC construction. This study employed the city of Nanjing as a case study and analyzed detailed information on how LULC changes over time of built-up areas. The surface area of the urban and built-up areas of Nanjing quadrupled from 11% in 1985 to 44% in 2015. In the same period, neither the entire city nor its subregions reached the VCR target of 80%. The proposed new methodology aims to support national, regional, and city governments to identify and prioritize where to invest and implement SCC measures more effectively in cities across China.


Introduction
With the implementation of the economic reform policy in 1978, the gross domestic product (GDP) of China has experienced an explosive increase along with urban expansion, especially in the eastern provinces. This rapid urban development has attracted a large number of the non-urban population, leading to massive urban migration. In turn, this has resulted in severe environmental damage, affecting human health and wellbeing of the Chinese urban population. The total amount of pollutants emitted by domestic and commercial activities, transportation, and industrial production has increased exponentially, exceeding the regional environmental carrying capacity [1][2][3][4][5]. A significant fraction of the emitted pollutants is transported with rainfall and runoff, having a severe chemical and ecological impact on the status of urban surface water bodies [3][4][5][6]. Moreover, urbanization leads to an increase of impervious area, which intensifies the hydrological cycle and increases the occurrences of extreme events such as droughts and floods [7][8][9].
Various concepts and solutions have been identified and implemented over the world, in order to reduce urban water pollution and to mitigate the adverse impacts of too much and/or too little rainfall, both in developed and developing countries. Some of them aim to expand the capacity of conventional gray infrastructure drainage systems, by using leaping weirs to intercept initial runoff, or by the construction of an underground storage facility to store large volumes of stormwater [10,11]. Green infrastructure is increasingly used to mitigate the impact of urban stormwater on both the quantity and quality of urban runoff. Approaches and strategies using green infrastructure to manage urban stormwater are referred to as best management practices (BMPs) in the United States and Canada [12,13], water-sensitive urban design (WSUD) in Australia and New Zealand [14], sustainable urban drainage systems (SuDS) in the United Kingdom and other EU countries [15][16][17], the Active Beautiful Clean Waters plan (ABC Water) in Singapore [18], and the Sponge City Construction (SCC) in China [19,20].
In terms of the design and assessment of Sponge City infrastructure, the volume capture ratio of annual rainfall (VCRa) is one of the compulsory indicators of the national SCC standard launched by the Ministry of Housing and Urban-Rural Development (MOHURD) in 2014 [20,21]. Examples of its application encompass Sponge City infrastructure projects covering an area ranging from a few to several hundred hectares [22]. The VCRa target value varies across China and is dependent on geographical factors such as climate, soil, and vegetation cover. Currently, a methodology is lacking to assess the spatial variation of VCRa at a city scale and changes therein over time, which hampers city planners to prioritize and select SCC interventions to comply with the SCC targets. This paper aims to present a method, based on land use and land cover (LULC) data derived from satellite remotely sensed images, to set a baseline of VCRa on city and district level required to identify, prioritize, design, and evaluate SCC projects by taking into account the VCRa in their pre-development stage.

Volume Capture Ratio of Annual Rainfall (VCRa)
The volume capture ratio of annual rainfall (VCRa) is defined as the fraction (expressed as a percentage) of the total annual rainfall that is not directly discharged (the cumulative annual control) but detained or retained by means of natural and artificial enhancement such as infiltration, storage, and evaporation [20].
The VCRa is a regional variant for calculating the water quality volume (WQv) in China [21]. Comparably, another similar practice by capturing the 85th to 95th percentile rainfall event has also been widely used [23][24][25][26][27]. The VCRa is based on the total captured stormwater volume of the survey site and does not consider the frequency of storm events [20,21,28]. Previous studies of VCRa mainly focused on three aspects: (i) regional application, (ii) case design application, and (iii) efficiency assessment of SCC measures or low impacted development (LID). With regards to its regional application, a map of China divided into regions with similar VCRa values was adopted in the SCC guide of the Ministry of Housing and Urban-Rural Development [20,29]. Li et al. [29] and Zhang et al. [21] proposed a design rainfall map corresponding to VCRa of 85% for the mainland of China. Although both maps are recommended for general application covering the whole country, the dataset used contains 186 rainfall stations [29,30], which is likely too limited to cover the more than 600 cities in China. Site-specific urban stormwater and flood models are commonly used methods for site design and assessment of VCRa. Liu et al. [31] developed an urban flood model (HIMS-URBAN) for the central area of Changde to design a Sponge City infrastructure project based on a short-duration design rainfall profile. Wang et al. [32] conducted a similar simulation using the Storm Water Management Model (SWMM) with gauged rainfall covering a period of three consecutive years. In both studies, a short-duration composite runoff coefficient was calculated to design the SCC project [31,32]. Randall et al. simulated four hypothetical scenarios including No LID (baseline), Low LID, Mid LID, and High LID measures for an urban catchment of 133 km 2 in Beijing to evaluate the feasibility of the VCRa target using SWMM [22]. Long-term (more than 30 years) continuous simulations were conducted on the basis of the conditions of the current landscape and stormwater infrastructure. On the basis of this study, the authors concluded that the computed VCRa served as a meaningful indicator to compare the effectiveness of the scenarios before and after the installation of SCC measures. In terms of efficiency assessment of SCC/LID measures, there is a wealth of publications focusing on the performance of LID measures and their design to meet water quantity and quality goals, such as from Hunt et al. [33][34][35], Cizek et al. [36], Braga et al. [37], and Drake et al. [38,39].
Various design manuals of SCC have been issued by the Chinese central government, and design targets have been set for different cities or regions [20,40]. For example, the VCRa of Nanjing for 2030 has been targeted at 80% [41,42]. A short-duration composite runoff coefficient is predominantly being used to assess the VCRa baseline on project site-level or district level [31,32,43], and the captured volume estimate is based upon urban stormwater management models [22,44]. Although methods to assess a short-duration composite runoff coefficient can be applied at the city level, a direct comparison with the design target of the VCRa is constrained due to the fact that the accuracy and precision of the underlying computer modeling of these methods are contested. This is partly due to the increasing complexity of urban stormwater simulation and the incompleteness or lack of high-resolution geospatial data affected by intense human activities [45]. The purpose of this study is to establish a method for a rapid global baseline assessment before the implementation of SCC measures.

Site Design
Nanjing, the capital city of Jiangsu province, is one of the megacities of the Yangtze River Delta Megalopolis in East China. Figure 1 depicts the study area with a radius of 20 km and a total surface area of 1256 km 2 , which covers a major part of the city of Nanjing. With a rapid economic growth, the GDP and population of Nanjing reached USD 191.1 billion and 8.245 million inhabitants, respectively, by the end of 2018 [46]. While economic growth has improved the livelihoods of the inhabitants, it has also led to urban congestion and a reduced capacity to manage extreme natural disasters. Since 2016, various green infrastructure projects and SCC technologies have been implemented in Nanjing to compensate against the environmental impacts of urbanization. VCRa targets of over 80% for new-built or refurbishment projects and over 70% for built areas have been adopted by the Nanjing city government. In addition, ambitious goals of "more than 20% of the built-up area have to comply with the requirements of the SCC targets by 2020" and "more than 80% will by 2030" were established [41,42]. Due to a lack of background information, it is difficult to assess the feasibility of these goals, which greatly reduces its guidance role for the SCC practices and measures chosen. The method designed in this paper is aimed to improve upon this gap. Compared to regular administrative division, a simple but effective division method has been used in Nanjing, whose administrative boundaries change frequently due to the rapid development of the economy and urbanization [48]. Four concentric circles with an increasing 5km step radius were drawn to analyze the relationship between landscape changes and distance from the center to outskirt, taking the Bronze Statue of Sun Yat-sen in Xinjiekou as the centroid to forming four concentric zones ( Figure 2). The surface area of the zones from the innermost to the outermost were 25π, 75π, 125π, and 175π square kilometers, respectively. From space, cities look like a kind of scratch on the planet surface, lacking uniformity or central symmetry. The concentric zones were divided by two crosslines through the centroid in order to capture the orientation of urban expansion. These two lines were set east-westward and north-southward to separate the areas of the concentric pie zones into quarters. This partitioning of urban area can also be applied to analyze other cities with a top-down planning and frequently changing administrative boundaries. Compared to regular administrative division, a simple but effective division method has been used in Nanjing, whose administrative boundaries change frequently due to the rapid development of the economy and urbanization [48]. Four concentric circles with an increasing 5km step radius were drawn to analyze the relationship between landscape changes and distance from the center to outskirt, taking the Bronze Statue of Sun Yat-sen in Xinjiekou as the centroid to forming four concentric zones ( Figure 2). The surface area of the zones from the innermost to the outermost were 25π, 75π, 125π, and 175π square kilometers, respectively. From space, cities look like a kind of scratch on the planet surface, lacking uniformity or central symmetry. The concentric zones were divided by two crosslines through the centroid in order to capture the orientation of urban expansion. These two lines were set east-westward and north-southward to separate the areas of the concentric pie zones into quarters. This partitioning of urban area can also be applied to analyze other cities with a top-down planning and frequently changing administrative boundaries.
Remote Sens. 2021, 13, x FOR PEER REVIEW 5 of 16 Figure 2. Partitioning of urban area as designed for Nanjing. The entire region was abbreviated as Z01, which covers the majority of the urban area. The concentric zones were distinguished with different colors (Z02-Z05). The concentric pie zones were annotated by text (Z06-Z21). Background from Google Earth Imagery [49]. The center was set at the Bronze Statue of Sun Yat-sen in Nanjing.

Data and Classification
The Landsat program has the longest continuous space-based record of the Earthʹs surface in existence [50,51]. On the basis of the service time and operational status of the Figure 2. Partitioning of urban area as designed for Nanjing. The entire region was abbreviated as Z01, which covers the majority of the urban area. The concentric zones were distinguished with different colors (Z02-Z05). The concentric pie zones were annotated by text (Z06-Z21). Background from Google Earth Imagery [49]. The center was set at the Bronze Statue of Sun Yat-sen in Nanjing.

Data and Classification
The Landsat program has the longest continuous space-based record of the Earth's surface in existence [50,51]. On the basis of the service time and operational status of the Landsat satellites, we chose a set of images to prepare the high-quality imagery for each year (Table 1) after the removal of the clouds and compositing operation. The Weka k-means clustering algorithm with 30 clusters was used to cluster the satellite image exporting preliminary classification [52][53][54]. Merging or reclustering and remapping the cluster numbers resulted in the 7 final land use maps. The classification was carried out using the Google Earth Engine platform [54].

State Vector, Transition Matrix, and Visualization
Mapping the spatial layout is traditionally based on a classification of land use. However, visualizing changes in land use remains a challenge. A method that combines the state vectors and transition matrix of land use and the Sankey chart was used in this study to analyze the process of urbanization in Nanjing. The transition matrix is a matrix to describe the changes in quantity from the former state to the latter state with the equations as follows: where P i = percentage of the total for land use type i at the former time, and ∑ m i=1 P i = 1; where P j = percentage of the total for land use type j at the later time, and ∑ n j=1 P j = 1; where P i j = percentage of the total from land use type i to type j; P i j also represents the percentage for which the land use type has not changed when i = j; the elements of the matrix match the following rules: The proportion of the part whose land use type has not changed can be calculated by the equation P unchanged = ∑ m i=1 P i i . Although it can have different numbers of classes for the different time points, it is recommended to keep them the same, i.e., m = n. In addition, the area of the specified land use type can be calculated by the area of the zone multiplied by P i or P j , and the area of change flows can be calculated by the area of the zone multiplied by P i j . The use of a Sankey chart to depict the state and dynamic process has been demonstrated in existing studies [55][56][57]. The original Sankey chart is good at visualizing the flows and relative visual proportions of components. Its disadvantage is that it has no axis, and thus the value or proportion cannot be read without the help of auxiliary text or annotations. In this paper, an improved method is presented, which combines the Sankey chart and a column chart to make it readable for both states and change flows. The modifications include the following: (1) a Cartesian coordinate plane is placed into the common Sankey chart where the x-axis represents time and the y-axis the relative contributions of the land use types; (2) for every specific time point, the number of nodes/columns is the number of land use types, and the height of them indicated the percentage of the total for the specific type; (3) the links between two time-verticals represent the transition matrix. This modified Sankey chart depicts the state components at certain time points as well as the change flows of the land use between those in one graph (Figure 3). For the calculation of the status vectors and the transition matrix, we used the Google Earth Engine platform [54], and for the visualizations of the modified Sankey charts, we deployed Observablehq (https://observablehq.com), an interactive JavaScript notebooks platform.

Precipitation Runoff Volume Control
The definition of volume capture ratio (VCR) is the ratio between the amount captured/retained and the amount of precipitation received at the site over a given period. The sum of VCRa and annual runoff coefficient (RC) equals 1 [58]. The upper limit of short-duration runoff coefficients for different land uses recommended by MOHURD [59] is listed in Table 2.

Precipitation Runoff Volume Control
The definition of volume capture ratio (VCR) is the ratio between the amount captured/retained and the amount of precipitation received at the site over a given period. The sum of VCRa and annual runoff coefficient (RC) equals 1 [58]. The upper limit of short-duration runoff coefficients for different land uses recommended by MOHURD [59] is listed in Table 2. For an urban zone composed of multiple land use/land cover types, a short-duration composite runoff coefficient can be determined by the area-weighted equation as follows: where RC c represents the short-duration composite area-weighted average runoff coefficient, and RC i , A i , and P i are the short-duration runoff coefficient, area, and percentage of the total area for the land use type i, respectively. Owing to the storage function of water bodies, the land use of water bodies is not considered part of the composed runoff coefficient. Moreover, it is widely recognized that the annual RC is smaller than the event-based short-duration RC due to the consideration of those rainfall events that do not produce any runoff [60]. As an inference, the short-duration RC can be considered as the upper limit of the annual RC, i.e., RC upper annual = RC event (5) where RC upper annual represents the upper limit of the annual runoff coefficient and RC event is the event-based composite runoff coefficient for each zone.
The lower limit of VCRa (VCRa − L) would be reasonable to depict the background value of the VCRa before the implementation of SCC measures, which can be calculated by the following equation: where VCRa − L, short for VCR lower annual , refers to the lower limit of the VCRa.

LULC Changed from 1985 to 2015 at the City Level
To avoid a massive pile-up of data and graphs, we highlighted three years in the results: 1985, 2000, and 2015 ( Figure 4). The states of land use in Z01 (which covers the majority of the urban area in Nanjing) changed considerably (Tables A1-A3, and Figure 5). The urban and built-up areas quadrupled from 11.09% of the total in 1985 to 44.07% in 2015, representing an absolute area size of 139 km 2 in 1985 to of 554 km 2 in 2015. The observed changes were mainly due to the rapid growth of the economy and urban population over the past few decades [61,62], a trend that has been seen everywhere in the territory of China. The area of water bodies increased by 20% in 2000 compared to 1985, which may be caused by the misclassification of the northwest shore of the Yangtze River. The water bodies in Z01 are mainly located at the Yangtze River and the Xuanwu Lake, which are protected areas for nature conservation to avoid its disappearance as a result of the rapid process of urbanization. The state values of mixed forests for the three years have not changed much, but the spatial distribution has changed greatly. This indicates that a transformation has occurred from a centralized forest or wood of the outskirt in 1985 to a built-up area with widely distributed greening trees in 2015 ( Figure 5). The croplands in Nanjing decreased irreversibly and the loss of its area was about 285 km 2 . The decrease of the croplands contributes to the majority of the increase of the urban and built-up lands. Grasslands made up the least area of land-use in Nanjing and its proportion remained relatively stable, while the spatial distribution changed considerably over the years. The reason may be that most of the grassland areas in the urban area are cultivated grass patches to serve the local needs of the citizens.  To distinguish the different LULC changes from the central district to the outskirts of Nanjing, we investigated the LULC changes of the concentric zones from Z02 to Z05 and the concentric pie zones from Z06 to Z21. The graphs of these results are shown in Figures A1 and A2.

5.2.
-L  To distinguish the different LULC changes from the central district to the outskirts of Nanjing, we investigated the LULC changes of the concentric zones from Z02 to Z05 and the concentric pie zones from Z06 to Z21. The graphs of these results are shown in Figures A1 and A2.

5.2.
-L Figure 5. The modified Sankey chart to visualize the state and change flows of the land use and land cover (LULC) classification for Z01 from 1985 to 2015 in Nanjing (colors and classes of the Sankey charts produced in this study are the same as in Figure 4).
To distinguish the different LULC changes from the central district to the outskirts of Nanjing, we investigated the LULC changes of the concentric zones from Z02 to Z05 and the concentric pie zones from Z06 to Z21. The graphs of these results are shown in Figures A1 and A2.

VCRα-L
VCRa − L, the lower limit of the VCRa, was designed as a baseline indicator to evaluate the pre-development stage of an SCC project. Figure 6 shows the spatial distribution of VCRa − L on the city-level from 1985 to 2015. Overall, 89% of the land area had a high VCRa − L of over 0.8 in 1985, and this proportion dropped to 54% in 2015. The numbers illustrate that more than half the area still meets the target of SCC in Z01. Nevertheless, the results of the composite VCRa − L for concentric pie zones showed a great difference. Overall, the value of VCRa − L increased with the distance away from the urban center. The smallest value always appeared at Z02 and the largest value always occurred at Z05 for these three time slices. On the other hand, VCRa − L decreased from 1985 to 2015. The values of VCRa − L in half of the zones were greater than 0.8 in 1985, and there were only three zones with a value below 0.5, which were Z06, Z14, and Z18. In 2000, two zones had an VCRa − L greater than 0.8; however, in 2015, none of them exceeded 0.8, and the eight zones were all less than 0.5. This indicates no subzone meeting the SCC target. The achievability varied for different subzones, i.e., Z18 is the most challenging subzone while Z21 is the easiest one, without consideration of the land area. High-density built-up in the center of Nanjing would be a tough project for SCC.
Remote Sens. 2021, 13, x FOR PEER REVIEW 10 of 16 -, the lower limit of the VCRa, was designed as a baseline indicator to evaluate the pre-development stage of an SCC project. Figure 6 shows the spatial distribution of -on the city-level from 1985 to 2015. Overall, 89% of the land area had a high -of over 0. 8 in 1985, and this proportion dropped to 54% in 2015. The numbers illustrate that more than half the area still meets the target of SCC in Z01. Nevertheless, the results of the composite -for concentric pie zones showed a great difference. Overall, the value of -increased with the distance away from the urban center. The smallest value always appeared at Z02 and the largest value always occurred at Z05 for these three time slices. On the other hand, -decreased from 1985 to 2015. The values of -in half of the zones were greater than 0.8 in 1985, and there were only three zones with a value below 0.5, which were Z06, Z14, and Z18. In 2000, two zones had an -greater than 0.8; however, in 2015, none of them exceeded 0.8, and the eight zones were all less than 0.5. This indicates no subzone meeting the SCC target. The achievability varied for different subzones, i.e., Z18 is the most challenging subzone while Z21 is the easiest one, without consideration of the land area. High-density built-up in the center of Nanjing would be a tough project for SCC. -for the concentric pie zones was annotated by text.

Conclusions
VCRa is a key performance indicator for the design and evaluation of SCC measures, but its assessment remains a challenge, particularly at the city level. A set of innovative methods using LULC data derived from remotely sensed Landsat images was developed and demonstrated in this research to conduct a baseline assessment (pre-development stage) of VCRa. For this purpose, a modified Sankey chart to visualize the states and change flows of LULC was developed. This method unifies mathematical descriptions and visualization of the state vectors and transition matrix to reveal the percentages and pathways of LULC changes. This novel and practical partitioning method has three advantages for its application in built-up regions: (i) it circumvents the variability of urban administrative boundaries, (ii) it enables a comparison between cities and districts within

Conclusions
VCRa is a key performance indicator for the design and evaluation of SCC measures, but its assessment remains a challenge, particularly at the city level. A set of innovative methods using LULC data derived from remotely sensed Landsat images was developed and demonstrated in this research to conduct a baseline assessment (pre-development stage) of VCRa. For this purpose, a modified Sankey chart to visualize the states and change flows of LULC was developed. This method unifies mathematical descriptions and visualization of the state vectors and transition matrix to reveal the percentages and pathways of LULC changes. This novel and practical partitioning method has three advantages for its application in built-up regions: (i) it circumvents the variability of urban administrative boundaries, (ii) it enables a comparison between cities and districts within a city, and (iii) it allows for a comparison of LULC between different historical time periods for cities with rapid growth. Moreover, a new ready-to-use indicator, VCRa − L, was proposed to evaluate the hydrological characteristics in the pre-development stage of builtup areas of a city and thus to set a baseline for the planning, design, and implementation of SCC projects.
A case study was conducted using the main built-up area of Nanjing (Z01) in China. These Sankey charts revealed the various LULC change pathways of urban landscape as follows: (i) significant changes of LULC occurred in Nanjing from 1985 to 2015, resulting from a rapid expansion of the built-up areas, a decrease of the croplands, and relatively limited changes of the surface areas of water bodies and green space; (ii) the built-up areas remained stable in the centrum zone (Z02), while the surface area of green space increased and of cropland almost disappeared, with the surface area of cropland decreasing significantly but still remaining important in the outskirt zone (Z05); (iii) preferential expansion of built-up areas occurred in the direction of southeast-northwest. From the perspective of VCRa − L, in 2015, more than half of the total city surface area still met the target of SCC, but none matched the target of 80% in the subzones.
The above case study further demonstrates that the modified Sankey chart method discloses LULC change pathways and patterns, which provide valuable historical information to understand the underlying causes of hydrological regime shifts required to identify realistic and effective SCC measures. The indicator VCRa − L can help policymakers and city planners to prioritize and select regions to implement SCC measures, to set realistic SCC targets, and to prioritize SCC projects. A shortcoming of the proposed method is that it does not allow an assessment at the level of a SCC project since its physical boundaries are set by functional or legal requirements instead of an arbitrarily chosen grid size.

Data Availability Statement:
The high-resolution LULC images for Nanjing at five-yearly intervals from 1985 to 2015 have been archived on Google Earth Engine Asset and Zenodo [63]. Moreover, the original data of the state vectors and transition matrixes of the LULC classification in Nanjing used to create the modified Sankey charts have been shared on Zenodo as well [64]. It would be easy to reuse for further analysis and will be distributed to those interested. . Figure A2. The modified Sankey charts to visualize the state and change flows of LULC for concentric pie zones (Z06-Z21) in Nanjing from 1985 to 2015.