Understanding spatial-temporal urban expansion pattern (1990–2009) using impervious surface data and landscape indexes: a case study in Guangzhou (China)

Abstract A new viewpoint for understanding the urban expansion using impervious surface information, which is obtained using remote sensing imagery is presented. The purpose of this study is to understand and describe the urban expansion pattern with the view of impervious surfaces instead of the conventional view of land use/land cover. Six years’ worth of impervious surface data (1990–2009) of Guangzhou are extracted via linear spectral unmixing analysis methods and spatial and temporal characteristics are discussed in detail. The area, density, and gravity centers changes of the impervious surfaces are analyzed to explain internal/external urban expansion. Meanwhile, five landscape indexes, such as patch density, edge density, mean patch size, area-weighted, and fragmentation index, are utilized to describe landscape changes of Guangzhou in past 20 years, which are influenced deeply by the impervious surface expansion. In order to detail landscape changes, two transects corresponding to the two urban expansion directions are designed and five landscape metrics in these two transects are reported. Conclusions can be drawn and shown as following: (1) temporally, the area of impervious surfaces increases from 12,998 to 59,911 ha from 1990 to 2009. The amount of impervious surface varies in different periods. The annual growth rates of impervious surface area during 1990–1995, 1995–1998, and 1998–2000 are 10.16%, 11.61%, and 10.78%, respectively; (2) annual growth rates decrease from 10.78% (1998–2000) to 5.67% (2000–2003). Nevertheless, from 2003–2009, the annual growth rate has a slight increase compared to a former period. The rate is 5.91%; (3) spatially, gravity centers of medium and high percentage impervious surfaces migrate slightly; and (4) according to the gradient analysis in the two transects, it can be observed that the high percentage of impervious surface increases gradually in new city districts (from west to east and from south to north).


Introduction
Cities represent one of Earth's fastest growing land-use types on a per-area basis, and over half of the planet's 7 billion humans now reside in cities. 1 It is very important to give a quantitative analysis in order to describe and forecast some problems that have occurred in the process of urbanization. Generally, remote sensing provided the most convenient means to monitor and quantify urban expansion by taking advantage of its wide area coverage and regular orbiting period and it had proven useful for mapping urban areas and obtaining data for the analysis of urban land cover change. [2][3] In China, urbanization (urban expansion) has been experienced during different periods during the last two decades and many projects or researches on urban expansion have been launched. [4][5][6][7] These researchers mostly placed their attention on land use/ land cover classification at the pixel level by using supervised or unsupervised classifications. These methods were limited by their precision of classification and image quality, such as spatial resolution, and they could not properly discriminate land cover classes in a heterogeneous context. 8 One of the recent major advances in urban expansion analysis was the vegetation-impervious surface-soil (V-I-S) model proposed by Ridd in 1995. In this model, impervious surfaces, which was defined as any surface which water could not infiltrate and which was primarily associated with transportation (streets, highways, parking lots, and sidewalks) and building rooftops, was the key indicator of city land cover/land use characteristics through the urbanization progress. 9,10 The V-I-S model was demonstrated to be an effective approach to cope with the mixed pixel problem and provided a guideline for decomposing low-resolution images of an urban surface and a link for these components to spectral signatures.
Meanwhile, many methods for extracting impervious surfaces were explored in recent years, such as decision tree method, supervised and unsupervised methods, maximum likelihood method, neural network, [11][12][13][14][15] regression tree model, 16,17 and spectral unmixing method. 18 Spectral unmixing analysis (SMA) had been widely utilized to derive the fractional contribution of representative urban land covers for its convenience and ability to derive physically interpretable information at a subpixel level. 19 The linear spectral mixture model was widely used for unmixing the mixed pixels. The spectral signatures in one pixel assumed a linear combination with their proportion as a weighting factor in the linear spectral mixture model. The linear spectral mixture model was an excellent approximation technique which was suitable for dealing with the spectral mixture problem, because it produced a set of maps that represented the abundance of each component.
The extraction of impervious fragmentation information was always conducted with the combination of high-albedo and low-albedo fraction images. 20 Although many researchers used SMA to extract impervious surface information, multitemporal studies for monitoring urban land cover changes were lacking. [21][22][23][24] Furthermore, there was no well-documented approach for analyzing the different internal and external urban expansion patterns through the urbanization progress.
Quantitative spatial analysis methods were needed in this research to explore the impervious surface sprawl development along the urban expansion directions. Gradient analysis, originally developed in the context of vegetation analysis, 25 had been used to investigate the effects of urbanization on plant distribution and ecosystem properties. 26 Luck and Wu 27 and Zhang et al. 28 used gradient analysis to study urban landscape patterns and the ecological consequences of the urbanization process by using a "split window" along their chosen transect. Kong and Nakagoshi 29 used gradient analysis combined with landscape metrics to characterize urban green space patterns quantitatively from 1989 to 2004 in Jinan (China) and gradient analysis was used to quantitatively understand the urban expansion of Guangzhou during the past 24 years based on the urban area extracted by the maximum likelihood method. 7 The primary objectives of this paper included: (1) examining and understanding the Guangzhou expansion progress and patterns using long time series impervious surface data (1990-2009); (2) analyzing detailed change information of landscape and gravity centers of impervious surface; and (3) explaining the internal and external urban expansion and sprawl patterns with the view of impervious surface under the subpixel level.

Study Area and Data Preprocessing
Guangzhou town is the core zone of Guangzhou municipality, which is chosen as the study area of this paper. Guangzhou municipality is the political, economic, educational, cultural, as well as scientific and technological center of the Guangdong province, South China. Guangzhou town covers approximately 3834.8 km 2 and is comprised of eight districts (Yuexiu, Dongshan, Haizhu, Liwan, Tianhe, Baiyun, Huangpu and Fangchun). Among these entirely city districts, Yuexiu, Dongshan, Haizhu, and Liwan are the traditional old districts, and the others are the new districts. It is bounded by the Pearl River to the east and south, Zhongshan city and Foshan city to the west, Dongguan city and Huizhou city to the east, and Qinyuan city to the north. There are three important highways, which are the Guang-Cong highway (from Guangzhou to Conghua), the Guang-Hui highway (from Guangzhou to Huizhou), and the Guang-Shen highway (from Guangzhou to Shenzhen).
projection system. The TM images of 1990, 1995, 1998, 2000, and 2003 are registered and rectified to the TM image of 2009 by the "image to image" method. After the geometric correction, the root-mean-square error (RMSE) is less than 0.5 pixels. At last, the vector data of Guangzhou administrative boundaries are applied to choose the study area subset from the whole image.

Masking Water Using Spectral Angle Mapper Method
Water's spectral characteristics were similar to the low-albedo urban surface; therefore, it was difficult to separate water from the low-albedo endmember, which would affect the unmixing results. A spectral angle mapper (SAM) was used to extract bodies of water and a mask was made to exclude pure water pixels from the Landsat images. SAM was a physically based spectral classification that used an n-D angle to match pixels to reference spectra. The algorithm determined the spectral similarity between two spectra by calculating the angle between the spectra and treating them as vectors in a space with dimensionality equal to the number of bands. Smaller angles represented closer matches to the reference spectrum. According to abundant experiments, which were design to explore the best value of the water spectral angle, the threshold identified was 0.1-0.08 rad. Only the pixels whose spectral angles fall into this range of radians were classified as water.
Water masking involved three steps: (1) selecting water endmembers; (2) utilizing SAM to extract the water and building a mask (spectral angle range is 0.1-0.08) via ENVI software; and (3) masking the water area.

Linear Spectral Mixture Analysis and Accuracy Assessment
Linear spectral mixture analysis (LSMA) was a linear unmixing method that calculated the abundance of the endmembers in a mixed pixel. To estimate impervious surfaces, the following four steps were generally involved: (1) reducing image noise via the minimum noise fraction (MNF) method, which produced a small number of components that contained most of the information; (2) calculating the pixel purity index (PPI) from the MNF component image. The PPI was the most representative interactive endmember extraction algorithm and was widely used due to its publicity and availability in the ITTVIS environment for visualizing images (ENVI) software originally developed by analytical imaging and geophysics; (3) selecting endmembers; and (4) unmixing the image based on the linear spectral mixing model.
The MNF rotation required two steps: (1) decorrelation and rescaling of the noise in the data based on an estimated noise covariance matrix, producing transformed data in which the noise had a unit variance and no band-to-band correlations, and (2) implementation of a standard principal component analysis of the noise-whitened data. The MNF rotation produced a two-part data set: one set was associated with large eigenvalues and coherent images, and another set was associated with near-unity eigenvalues and noise-dominated images. In MNF rotation, noise was separated from the data and spectral processing results were improved. 19 Using the majority of the information from the original image, the first three components were used to select endmembers. The selection of endmembers was complemented with the use of the PPI function by projecting three-dimensional (3-D) scatter plots of MNF fraction components 1, 2, and 3 onto a 2-D space and marking the extreme pixels in each projection. Pixels with a high digital numbers (DN) value indicated the spatial location of spectral endmembers. 30 The mathematical model of LSMA could be expressed as where i ¼ 1; : : : ; m (number of spectral bands); k ¼ 1; : : : ; n (number of endmembers); R i was the spectral reflectance of band i which contained one or more endmembers; f k was the proportion of endmember k within the pixel; R ik was the known spectral reflectance of endmember k within the pixel on band i; and ER i was the error for band i. A constrained least-squares solution was used in this research, assuming that the following two conditions were simultaneously satisfied: Endmember selection was a key step in the LSMA approach. Endmembers could be collected through laboratory or field measurements, or from images. 31 According to the V-I-S model, three basic endmembers were chosen: vegetation, impervious surface, and soil; the vegetation was the sum of forest, farmland, and grassland fractions. In addition, the impervious surface's extraction had always been conducted by the low-albedo and high-albedo fractions. Because the lowalbedo fraction image was more complex than other fraction images, the direct extraction of impervious surface image from the high albedo and low albedo fractions was obviously unsuitable. Wu and Murray (2003) built a model to extract the impervious surface and this idea was accepted widely. 20 The model is shown as follows: where R imp;b ,R low;b , and R high;b were the reflectivities of impervious surface, low albedo, and high albedo of Band b. f low and f high were the densities of low albedo and high albedo. e b was the value of the remaining data. Equation (4) must meet the needs of the following equation: Assessment of impervious surface images was important for the validation of classification results. Some researchers had discussed the accuracy assessment approaches, 32 including overall accuracy (OA), producer's and user's accuracy (PA and UA) and kappa coefficient, and so on. Much research used RMSE to assess the results in the accuracy assessment procedure. 15,[33][34][35] RMSE and systematic error (SE) were utilized in this research to evaluate the accuracy of impervious surface estimation. RMSE and SE could be expressed as follows: where X ͡ i was sample i's density of the impervious surface estimated from TM and ETMþ data, X i was sample i's corresponding value of ground reference images, and N was the number of the sample.

Gradient Analysis
In this research, two transects ( Fig. 1) were designed to measure the impervious surface sprawl. The principle of design was listed in following steps: (1) for the whole study area, the impervious surface was extracted first using the LSMA method, (2) the impervious surface value of 2009 were analyzed using Statistics tools of ENVI software, and the pixels, which had highest value, were selected. (3) if the number of pixels with the highest value was one, then this pixel was considered as the origin point. If the number of pixels with highest value was greater than 1, the 8 neighborhood pixels were computed and the mean value of these 9 pixels was obtained. The pixel with the highest mean value (9 pixels) was selected as the origin point; (4) by such an analogy, the selection process would stop until only one pixel was selected; (5) after the origin point was selected, 12 directions every 30 deg were checked based on this origin point. These two directions had the longest distance with the largest sum of impervious surface (Fig. 1); and (6) along with these two directions, two transects were designed. These two transects were included with a square that was called "Cell." It turned out that these two transects extend from Guangzhou New Baiyun International Airport in Baiyun to southern Fangcun, and west from Liwan to eastern Huangpu, and cut across the entire study area. In this paper, each transect was divided into 36 cells (each cell covered 1 × 1 km approximately, i.e., each cell was made up of 33 × 33 pixels). From 1990 to 2009, the expansion distances were 21 km (south to north) and 25 km (west to east), respectively, and the annual expansion distance was about 1 km (i.e., Guangzhou expands 1 km outward every year). Therefore, the scale of 1 km conformed to the real expansion state of Guangzhou. Therefore, the 1 km scale was made as a cell to discuss the spatial-temporal change of landscape metrics.

Landscape Indexes
To quantify the spatial pattern of exterior urbanization along the two transects, five landscape indexes in each cell of these two transects were calculated using patch analysis tools of GIS software. These five landscape indices were patch density (PD), edge density (ED), mean patch size (MPS), area-weighted (AWMSI), and fragmentation index (FI). In this paper, the landscape meant "cell," which was designed in Part 3.3 above, and the "patch" was the set of "pixels" in each cell. These pixels were selected according to given value of impervious surface area.
PD was a limited but fundamental aspect of landscape structure. PD had the same basic utility as the number of patches as an index, except that it expressed the number of patches on a per unit area basis that facilitated comparisons among landscapes of varying sizes where n i meant the total area of impervious surface and A meant the total area of each cell. ED at the class level had the utility, limitation, and standardized edge to a per unit area basis that facilitated comparisons among landscapes of varying sizes and equaled the sum of the lengths of all edge segments involving the corresponding patch type where E meant the total length of the edge in pixels and A meant the total area of each cell. MPS was another class and landscape index based on the number of patches. A landscape type with a smaller MPS than another patch type might be considered more fragmented. Thus, MPS could serve as a habitat fragmentation index where a ij meant the area of impervious surface and n i meant the number of pixels in one cell. AWMSI was a useful form of the shape metric. With AWMSI, the shape of each patch was weighted by its area relative to the area of the corresponding size. Specifically, larger patches were weighted more heavily than smaller patches in calculating the average patch shape for the landscape.
where P ij meant the perimeter of the pixels and a ij meant the area of the pixels.
FI characterized the degree of fragmentation of the landscape that reflects the complexity of the spatial structure of the landscape where n i meant the number of pixels in a cell and A i meant the total area of a cell.

Area and Rate of Impervious Surface
According to the V-I-S model, four endmembers are selected, which are high albedo, low albedo, vegetable, and soil. Using the LSMA method, four feature maps of high albedo, low albedo, vegetable, and soil are generated. Impervious surface maps are yielded by the method that is discussed above (low-albedo image þ high-albedo image). The impervious surface area of each year is extracted (Fig. 2). The RSMEs of the extraction results are 0.16, 0.20, 0.21, 0.31, and 0.19, respectively. Areas, increased area, and growth rate (1990-2009) of impervious surfaces are listed in Tables 1 and 2.
As shown in Tables 1 and 2, the impervious surface area increases from 12,998 ha in 1990 to 59911 ha in 2009. This implies that the impervious surface has undergone a dramatic development in the past two decades. Although the annual growth rate of impervious surfaces became lower than in the past after 1998 (10.78%), the impervious surface density of Guangzhou still kept growing. The most important factor to explain this phenomenon is that in 1990, Guangzhou was in the initial stage of urbanization. In this stage, there were few impervious surfaces (built-up areas), and a lot of vegetation, farmland, soil, and bare land were the main land cover components. With the development of the economy, Guangzhou went into the advanced stage of urbanization. A large amount of vegetation, farmland, soil, and bare land changed to impervious surface area quickly and inefficiently. With the passing of time, there were no more spaces to support the extensive urban expansion. In order to satisfy the demand of urban development, the density of impervious surface kept growing in the inner city, and the utilization efficiency of land was higher than before.

Spatial Migration of Gravity Center (Medium and High Percentage Impervious Surface)
The impervious surface is classified into three categories based on its percentage: (1) low percentage impervious surface, the value of impervious surface is 0%-30%; (2) medium percentage impervious surface, the value of impervious surface is 30%-70%; (3) high percentage impervious surface, the value of impervious surface is higher than 70%. According to Fig. 3, it is obvious that in 1990, 1995, and 1998, the medium percentage impervious surface increased year after year. In 1995, the medium percentage impervious surface increased sharply by almost 15,000 ha, which is mainly due to spatial sprawling of the impervious surface area in the first stage of urbanization. Meanwhile, high percentage impervious surface area appeared and increased gradually in this period. The reason is the rapid development of the central district of Guangzhou after the reform and opening-up policy. In this period, the strong demand for land and environment/transport facilities arose dramatically, thus the city expanded outward quickly. A great decrease in medium percentage impervious surface area (almost 9000 ha are lost) came out in the period from 1998 to 2000 and a steady growth occurred in the period from 2000 to 2009. The high percentage impervious surface increased sharply from 1998 to 2000, that is to say, parts of medium percentage impervious surface were developing into high percentage impervious surface area step by step. Afterward, the high percentage impervious surface area experienced a high growth rate and finally surpassed the medium percentage impervious surface area in 2009. This implies that, in the past two decades, the external expansion of impervious surface area dominated the early urbanization. When the impervious surface area expanded to a certain degree, there was no more space to support its expansion; therefore, impervious surface area changed into an internal expansion pattern; in this stage, the percentage of impervious surface area became higher than that in the first stage.
Furthermore, the spatial-temporal evolution pattern can be understood very clearly on the spatial aspect through gravity center migration. The gravity center in spatial analysis is borrowed from Newton's famous law of gravitation. The gravity center is identified as the point where the equilibrium of the body can be reached. In spatial analysis, a gravity center of a region is the Fig. 2 The impervious surface fractal images of 6 years (1990, 1995, 1998, 2000, 2003, and 2009 location where the maximum accessibility can be reached. The maximum accessibility means that the gravity center has the total minimal distance to all other locations. Finding the gravity center can optimize the process of site selection, which can help one to find the relatively best site for urban development. In this paper, the gravity center is used for finding the best site, which can represent the impervious surface of the whole study area. In ArcGIS, there is a "Spatial Statistics" geo-processing toolbox with a tool called "Mean Center," which can identify the geographic center or the center of concentration for a set of features. The gravity centers of medium and high percentage impervious surface in each year are calculated (Fig. 4).
The gravity center of medium percentage impervious surface area moves from the west of Baiyun Mountain toward northeastern of Guangzhou (a-b-c-d-e-f). The movement distance of the gravity center is almost 10.67 km. The annual movement rates in different periods are calculated to reflect the different paces of medium and high percentage impervious surface along the urban expansion process, respectively. The mean rates of medium percentage impervious surface gravity center in different periods are 0.   The density of impervious surface in the internal city has increased since 1995 and has gradually become saturated. Spatially, the gravity center of high impervious surface appears in Tianhe district in 1990. After that, it moves to northwestern (Baiyun district) in 1995, back to the Tianhe district in 1998 and back to the Yuexiu district in 2000 because of the commercial/industrial developing of Tianhe and Yuexiu. It moves back to Baiyun district once again in 2003 because of the construction of New Baiyun International Airport and the facilities along the airport roads. In 2009, the gravity center moves to northeast of Guangzhou, which means much more high percentage impervious surface is appearing in the junction of the Tianthe and Baiyun districts.

Gradient Analysis with Landscape Indexes
It is found that the important characterization of urbanization is the dramatic increase of medium and high percentage impervious surface. In this section, impervious surface area percentage higher than 30% is chosen to calculate the landscape indices in each cell (1 × 1 km size) to investigate the urban expansion pattern. In this paper, five landscape indices are calculated using patch analysis of GIS software. Results are shown in Fig. 5.
As shown in Fig. 5 chart (a), PD values increase gradually, with the distance from the original point becoming further year by year. The PD reaches the first peak values in cells 7-9, which are located in the Tianhe district, and the second peak value in cells 23-29, which are located in the Huangpu district in 2009. In 1990, the highest value of PD appeared in cells 5-6, which are located in the core of the city (inner city). However, in 2009, the highest value of PD appears in cells 25-26 (external city). The value of PD in 2009 is almost double that in 1990. In the other transect, chart (b), the same situation also occurs as well. In 1990, PD reaches its first peak value in cells 5-6 (inner city) and second peak in cells 18-19 (external city). Based on the PD changing in these two transects, it is clear that Guangzhou has had an external expansion process in past 20 years.
Charts (c) and (d) show that the MPS increased year after year. The same situations can be found from these two transects in the 6 years that the peak of MPS is located in cells 0-9 (inner city) and the lower values of MPS are located in cells 17-30 (external city). The phenomenon declares that, in the old city, vast stretches of impervious surface have been forming with urban development. According to charts (e) and (f), in two directions the ED, which reflected the degree of fragmentation of the impervious surface, has a maximum value in 1990. However, ED has a minimum value in 2009 compared to other years. Therefore, the conclusion is drawn as follows: impervious surfaces are most fragmental in 1990 and most continuous in 2009 in these two directions.
It is obvious that in charts (g) and (h), AWMSI increases gradually while the distance from center and ranks decline in both directions in 1990. At the same time, AWMSI value varies irregularly year after year. In charts (i) and (j), FI ranks extremely high in 1990 and increases with the distance from the city in every year. In 1990, FI keeps a steady value around cells 5-15 in the core city. In chart (i), FI reaches the peak value at cells 25-26 until 1998, which is located in New Baiyun International Airport and FI decreases gradually from 2003 because of the establishment of New Baiyun International Airport (the establishing time is from 1998 to 2003). In chart (j), FI reaches its peak value for each year at cells 14-18, which are located in the new city. Since 1995, the FI increases year by year which indicates the medium-high density impervious surface begins in two districts. After 2003, the FI decreases because the urban development is stable.
Near the center of city, along the northwest-southeast direction from 3 km (left) to 3 km (right) and the southwest-northeast direction [3 km (left) to 3 km (right)], values for PD show a visible decrease between 1990 and 2009. Meanwhile, values for MPS show a different temporal pattern. It has a progressive increase in the period from 1990-2009. It can be inferred that the impervious surface around the city center developed from several small sized impervious surface areas into fewer larger sized impervious surface areas. Values of ED and FI, near the city center, along the northwest-southeast direction from 2 km (left) to 2 km (right) and the southwest-northeast direction [3km (left) to 3km (right)], have an apparent decrease during the period from 1990-2009. This indicates that the impervious surface area distributed around the city center evolved from very fragmental to more continuous. Considering the values of PD, MPS, ED, and FI together, these patterns suggest that around the city center, impervious surface area developed from sparse little sized patches to a dense distribution with a larger area. This is due to the dispersed urban sprawl that characterized these areas starting in 1990, and the "fillingin" process that subsequently occurred, resulting in a higher residential density distribution.
Near the fringe of the study area, in the northwest from 4 km (left) to 3 km (right) and southeast from 3-31 km (right), the peak values of PD and MPS in each of the 6 years reflect that impervious surface area was gradually sprawling from the internal city center to external city districts after 1990. Subsequently, some new centers appear with a higher density distribution of residential and industrial buildings. It is inferred from this that impervious surface area appeared irregularly and sparsely near the edge of the city along the northwest-southeast direction and indicates that the incipient urban external sprawl is fast and unplanned. In addition, in the southwest, from 9 km (left) to 3 km (left) and northeast [between 3 and 27 km (right)], the same pattern with a northwest-southeast direction has been found, which exactly reflects the noncenter urban spatial sprawl pattern.

Conclusions
In this study, the spatial-temporal patterns of urbanization in Guangzhou have been investigated by using a subpixel linear spectral unmixing method in the periods from 1990-2009. Gradient analysis is finished based on two mutually perpendicular transects along the urban expansion directions, including 72 cells combined with landscape metrics. Impervious surface area is obtained and analyzed for 6 years. Landscape indexes are calculated based on the distance from the original point for medium and high percentage impervious surface to analyze the urbanization spatial-temporal pattern along the urban sprawl directions.
The results indicate that Guangzhou has experienced a rapid expansion in the past two decades (1990)(1991)(1992)(1993)(1994)(1995)(1996)(1997)(1998)(1999)(2000)(2001)(2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009)). The impervious surface area increased from 12,998 ha in 1990 to 59,911 ha in 2009 spatially. In the whole study area, the gravity center of medium percentage surface gradually migrated from west of Baiyun Mountain toward the northeastern direction. The reason is that other types of land cover located in northeastern Guangzhou has been converted into impervious surface area. The gravity center of high-density impervious surface migration leapfrogged, for the reason that human activities and policies more severely affect the high-density impervious surface development. Along the urban expansion axes, medium and high percentage impervious surface appears gradually in the new city districts (noncenter area), from a sparse and irregular pattern to a less fragmental and dense distribution, which indicates that dispersed urban sprawl characterizes these areas in 1990 and the "filling-in" process occurs subsequently, resulting in a higher residential density distribution.
The expansion pattern of impervious surface is an important basic ingredient for analyzing urbanization. Linking the quantified and spatially explicit impervious surface sprawl pattern can help to clarify the relationship between the urbanization process and impervious surface expansion. Finally, developing an understanding of the dynamic spatial patterns of impervious surface area can improve our ability to assess and create future planning scenarios by combining appropriate methods, such as gradient analysis. The ability to relate pattern to process is essential for understanding urbanization, and this cannot be done effectively without gravity center and landscape metrics analysis.