Abstract

In this work, we present a processing chain for landscape pattern and ecological security status assessment and prediction based on cellular automata Markov (CA-Markov) and pressure status response pattern (PSRP) models using remotely sensed data (RSD) captured in 1986, 1996, 2006, 2016, and RSD simulated in 2026 over Zhengzhou city, Henan province, China. Three major findings can be withdrawn through the experiments. First, there is a significant changing of landscape type area, especially for building land. The area of building land is up to more than 5%, from 1986 to 2016. Secondly, the heterogeneity of landscape is increasing, and the diversity of landscape is becoming more and more diversifying and complex. Third, the changing trend of ecological security of Zhengzhou city shaped as decreasing and increasing gradually during the last 40 years. While the ecological security status, nowadays, appeared to a good trend by contrast of the previous stages. The predicted results with CA-Markov model show that the level of ecological security is still in moderate and has a trend of moving toward to the center in 2026.

1. Introduction

The change of land use and land cover (LULC) is the key point for urban sprawl research and is becoming an advantage and hot pot for research of global environmental change and sustainable development [1]. The changes of landscape, as a complex structure of land use, at different scales, shapes, and quantities, express inherence mechanism of land use and land cover change (LULCC) [2]. The research of interaction mechanism between landscape pattern and LUCC has become an essential scientific topic in earth surface environmental change [3]. Ecological security assessment (ESA) is the foundation of regional ecological environment governance and prevention, because the comprehensively ecological risk and vulnerability of natural, economic, and social levels are considered in the changing procedure of landscape pattern [4]. And the integrity of ecosystem and health of overall ecosystem status can be reflected by ESA.

At present, more and more attentions are paid to the changing trend, vulnerability, sensitivity, and heterogeneity analysis methods of landscape pattern, such as landscape pattern index, the particle size effect, and spatial statistics algorithm [2, 57]. The wildly accepted ecological security assessment system include the “ecological footprint” concept and model, the concept of pressure-state-response (PSR) framework model, and the concept of nature-economy-society framework model [8, 9].

The motivation of this paper is to develop a processing chain for landscape and ESA monitoring, analysis, and prediction with mathematical and geography models. The objective of this work was to study LULC change statues during 40 years and to research relationships among LULC and landscape pattern and ecological security.

The paper is organized as follows: a brief introduction about research area and preprocessing of datasets are given in Section 2. The proposed chain based on remotely sensed observation techniques and CA-Markov model, PSRP model, and analytic hierarchy process (AHP) method are applied to the test area in Section 3. The results and analysis are shown in Section 4, and the conclusions are drawn in Section 5.

2. Study Area and Datasets

2.1. Study Area

Zhengzhou city, capital of Henan province, China, which is situated at the transitional zone among piedmont, torrential plain, and the Yellow River alluvial plain, is selected as the study area. The geographic coordinates lie in 112°42–114°14E, 34°15–34°57’N, shown in Figure 1. The city is bounded by Nanyang Road. West part of the road is the piedmont inclined plain with the height 100–130 m and surface lithology of Pleistocene series loess powder on the soil. Moreover, east part is the Yellow River alluvial plain with low height 85–100 m. This test area has climatic conditions of continental monsoon climate, with a total area of 7446 square kilometers and 11 million resident population. The annual average temperature is 14.5°C. The general terrain trend is tilt from southwest to northeast. Nowadays, the test area is the center of politics, economy, science and technology, culture, finance, transportation, and information in Henan province and undergoing Chinese urban agglomeration, and economic development and urban expansion are accelerating.

2.2. Datasets and Preprocessing

The datasets collected in this research include multitemporal remote sensing (RS) datasets, LULC product data, simulated LULC data in future, vector data, and demographic data. The selected multitemporal remote sensing (RS) datasets acquired for four different times include Landsat TM/ETM+ (captured on 12 August 1986, 15 April 1996, 9 April 2006, row in 36, and path in 121), Landsat OLI (captured on 23 April 2016, row in 36, and path in 121), GlobeLand30 product obtained in 2010, simulated LULC map in 2026, the division vector data of administrative area with a scale of 1 : 30000, and population statistic data (captured in 1986, 1996, 2006, and 2016, respectively). More detailed description about datasets collected can be found in Table 1. The data availability statement: all datasets used in this paper including the primary remotely sensed data and boundary vector data, as well as process results data can be obtained from https://pan.baidu.com/s/1rThlB5U2RQExrKAlM26Zig. Researchers who are interested in this topic can download the data to replicate the analysis or conduct secondary analysis.

The preprocessing includes geometry correction, remotely sensed data resizing, and LULC types obtaining. In detail, (1) the geometry correction is applied to different satellite images based on the topography map over test area using the polynomial transformation equations. The degree is set to one, and the nearest neighbour resampling method is selected in the polynomial method, which is defaulted by ENVI® software. In this process, firstly, seventy-four ground control points, seventy-seven ground control points, and seventy-eight ground control points are selected between remotely sensed data 1986 and 2016, 1996 and 2016, and 2006 and 2016, respectively. Then remotely sensed data obtained on 1986, 1996, and 2006 are geocoding to the same geometry system based on remotely sensed data obtained on 2016. The root mean square error (RMSE) value for different scenes is less than 0.8 pixel (0.67, 0.71, and 0.56 for 1986, 1996, and 2006, resp.); (2) the 30-meter resolution remotely sensed data over test area in 1986, 1996, 2006, and 2016 are obtained with data resize method based on administrative boundary vector file. In this process, the spatial subset of remotely sensed data over research area was obtained through data resize method based on vector data. Firstly, the vector data was imported and overlapped to the whole remotely sensed data. Then the resize model with spatial subset in ENVI packages was performant.

3. Methodology

3.1. LULC Map Classification

LULC status are categorized into five classes using the support vector machine (SVM) classifier [10] according to the Classification Criteria for Land Use Status/GB-T21010–2015 and GlobeLand30 standard products [11]. The GlobelLand30 is an open-access 30 m resolution global land cover data product which comprises ten land cover types, including water bodies, wetlands, artificial surfaces, cultivated lands, forests, shrub lands, grasslands, and barren lands [12]. Considering the situation and comparing with GlobelLand30 product of the selected research area, the selected five LULC classes are agricultural land, vegetation, water body, bare land, and building land. The selected SVM classifier, which uses a various of kernel function to construct nonlinear decision boundaries in the original feature space and transfers the sample data into linear boundaries of a higher dimension, has been successfully applied for regression, classification and has been demonstrated as one of the most powerful and efficient tool for remotely sensed data classification [13]. Therefore, in this research, the SVM classifier which derived from statistical learning theory was selected to obtain LULC status. The classification scheme proposed in this research are as follows. Firstly, about 100 points of training samples for each class are selected by expert in remote sensing field; secondly, the radial basis function (RBF) kernel, where parameters are determined automatically to minimize the upper limit of the expected test error, and according to the results from previous research, it works well in most cases, was selected in this work. The mathematical representation of RBF is , where is the bias term in kernel function for polynomial and sigmoid kernels; third, inputting feature vectors, such as band 1 to 5 and band 7 for Landsat TM/ETM+ data; and lastly, applying SVM classifier to feature combinations obtained in the third step and obtaining LULC results. Additionally, LULC map of 2026 was simulated based on LULC map of 2006 and 2016 using CA-Markov model. The thematic maps produced using SVM classifier in different years were accuracy evaluated for further procedure, such as change trend analysis, LULC trend prediction, landscape pattern and ecological security assessment, and so on. The commonly used accuracy evaluation indexes including overall accuracy (OA), user’s accuracy (UA), producer’s accuracy (PA), and kappa coefficients (KC) can be performed based on matrix statistics [31modelling]. In this research, the overall and comprehensive indicator include OA and KC, which can be calculated using expressions and , respectively, are selected. Where is number of pixels correctly classified, is the total number of pixels in the confusion matrix, is the number of rows, and and are the column and row total, respectively [14]. There are total of 15,538, 13,829, 15,991 and 13,937 pixels selected for training in remotely sensed data captured on 12 August 1986, 15 April 1996, 9 April 2006, and 23 April 2016, respectively. The number of selected test pixels (ground truth points) are 5643 points. The classification map results are shown in Figure 2. The classification accuracy shown in Table 2 demonstrates that LULC maps obtained with our method are good enough for further processing.

3.2. SDDLUC and IDDLUC Model

The changing speed and changing amplitude are analyzed by calculating landscape type transfer matrixes in order that changing trend can be characterized from spatial and temporal aspects. The variation of amplitude is expressed in the area change of landscape, while the changing speed is embodied in the single dynamic degree of LULC and integrated dynamic degree of land use change [15, 16]. The single dynamic degree of land use change (SDDLUC) describes the changing degree rate of a specific land use type of the area. SDDLUC can be described in the following equation [17]: where is the SDDLUC index, and are areas of a specific land use type at the beginning and ending time, respectively, during the study period. is the length of the research period.

The integrated dynamic degree of land use change (IDDLUC) describe overall changing degree of the study area with a certain period. IDDLUC can be described in the following equation [18]: where is area of the land cover type in initial time, while is transfer area amount of the land cover type to other land cover types during the research period.

3.3. Simulation with CA-Markov Model

The CA-Markov is a combination model of Markov model which describes the probability of land cover changes between start and end period by developing a transition probability matrix between them and CA model which allows the transition probabilities of one pixel to be function of adjacent pixel [14]. The CA-Markov model has a capability to simulate and to predict LULC changing status in future with remotely sensed data [19]. It is helpful to the research area in establishing and optimizing urban development decisions at different spatial-temporal dimensions [14, 1921], when simulating and prediction urban landscape changes comprehensively in quantity and spatial aspects. As given in [19] by Singh et al., the Markov chain process can be described as , where means Markov chain process for a particular time (t),, , and means previous time, current time, and end time in future, respectively.

In order to realize this objective, strategies can be divided into several procedures as follows: (1) setting LULC classification results obtained using classification method as initial variable; (2) calculating transition probability matrix and condition probability matrix of LULC types as conversion rules with GIS analysis method; (3) analyzing accuracy of simulating and predicting results based on the contrast of classification maps and predicted maps; and (4) changing trend analysis and ecological security assessment on the basis of classification results and simulation results. In this particular experiment, the four scenes LULC maps obtained using SVM classifier from remotely sensed data captured in 1986, 1996, 2006, and 2016 were set as initial variable and used to develop transition probability matrixes for LULC types between adjacent year; and then 10 times for CA iteration and 55 standard filter was selected to simulate LULC maps in 2006, 2016, and 2026, respectively. And simulation accuracy was calculated through contrasting and comparing the simulated LULC maps in 2006 and 2016 with classification maps obtained in the correspondent year. Then the changing trend and ecological security status between the periods of 1986 to 2026 are analyzed and assessed.

3.4. Landscape Pattern Index Analysis

The landscape pattern indexes are important indicators to quantitative describe and monitor the transformation characteristics of landscape structure. While it is a challenge to choose appropriate landscape indexes that would reflect the landscape properties and is not necessary to calculate all landscape indexes for a special landscape, since some of these indexes are strongly correlated. The choice of landscape pattern indexes can be based on both quantitative analysis and knowledge of the study area [19, 22, 23]. One of previous study suggests the combining use of Spearman’s correlation coefficient (SCC) and principal component analysis (PCA) method to reduce redundancies [22]. The method proposed in [22] used SCC to reduce high correlation between two landscape indexes by checking SCC whether greater than 0.9 or less than −0.9 and used PCA [14, 2224] to reduce the number of landscape indexes to a manageable level. In this research, the selection of appropriate landscape indexes was based on experience of previous research as well as the knowledge of the study area when the spatial distribution status and data quality of research area were considered. The slight correlation landscape pattern indexes at both patch level and landscape level are selected to analyze long-term series of landscape pattern change of research area from 1986 to 2026 [25]. The number of patches (NP), patch density (PD), largest patch index (LPI), edge density (ED), mean shape index (MSI), fractal dimension index (FDI), aggregation index (AI), and interspersion juxtaposition index (IJI) were selected when landscape pattern was analyzed in patch level. While number of patches (NP), patch density (PD), largest patch index (LPI), fractal dimension index (FDI), contagion index (CONTAG), Shannon’s evenness index (SHDI), Shannon’s diversity index (SHEI), fragmentation index (FI) were selected when landscape pattern was analyzed in landscape level. The integrated application of multiple landscape indexes can avoid the limitation in some extent when individual landscape index was selected and can reduce the adverse effect on result analysis caused by remotely sensed data resolution and human factors.

3.5. Establishment of ESA System

The ecological security assessment system established in this research is based on the concept of pressure-state-response (PSR) put forward by The UN Organization for Economic Cooperation and Development (OECD). The establishment of PSRP model is according to the principle of dynamism, hierarchy, and feasibility [9, 2628], and according to the stability of the ecosystem structure, the resilience of the landscape, the service function of ecosystem, and the changing factors of landscape pattern. The established so-called pressure-state-response-pattern (PSRP) model system included target layer, criterion layer, index layer, and 12 single indexes [29, 30].

While there is no comparability between the various indexes because of different units of different indexes, for the sake of quantification and computation, the primitive value of indexes is standardized as between 0 to 10. The standardization model for positive indicators (the larger the values, the safer the ecology) can be described in the following equation:

The standardization model for negative indicators (the smaller the values, the safer the ecology) can be described in the following equation: where is standardization value, is primitive value of indexes, is the maximum value of a set of indicators, and is the minimum value of a set of indicators.

Different indexes have significant different contribution for ecological security assessment. To overcome the diversity and to guarantee the relative importance of indexes in the ecological security assessment system, entropy method and AHP method were applied comprehensively to determine weight of indicators (results are shown in Table 3). Entropy method was used to determine the weight coefficient of indexes by using the diversity of the difference of indexes information, which can avoid the randomness caused by subjective decision. At the meantime, as a complement, AHP method was employed to verify the results of weight determined and to overcome the limitation of absolute objectivity in entropy method [31].

With the help of GIS spatial analysis mechanism, the matrix of transition probability and matrix of condition probability were calculated as conversion rules in evaluation level of ecological security assessment. The status of cellular automaton model is defined as the class level of pixels for ecological security assessment. Taking the data in 2006 and 2016 as an instance, ecological security assessment results in 2006 and 2016 were setting as initial status, with standard 55 filter and 10 times iteration parameters, the ecological security trend can be predicted [21].

To make sure the equitable using various indexes in assessment system and establish the importance of different indexes, a linear weighted synthetic evaluation model (LWSEM) was built. The LWSEM can be described in the following equation: where is comprehensive evaluation value of ecological security; is weight of index; is standardized value of index. Domestic and foreign methods of classifying level of ecological security were referenced to determine criteria of ecological security assessment class [26]. The established criteria of ecological security assessment are shown in Table 4.

4. Results and Analysis

4.1. Changing Speed and Amplitude

To analysis landscape changing character, LULC changing speed and amplitude between 1986 and 1996, 1996 and 2006, and 2006 and 2016 need to be calculated in advance. The results were shown as Tables 57.

The results of LULC changing statistics shown that the transformation relationships of various landscape types were complex, and the variation range of the single dynamic degree was wide in this research area. While the integrated dynamic degree during three periods (1986–1996, 1996–2006, and 2006–2016) was 21.65%, 20.17%, and 19.21%, respectively. The integrated dynamic degree indicated that there was a dramatically landscape type change in this period. The characteristics of landscape changing were discussed as follows: (1)There was a downward-downward-growth changing trend of water body. There was a sustained downward trend of water body from 1986 to 1996 and from 1996 to 2006. The changing amplitude of water was great, and especially during the period between 2006 and 2016 was the biggest among all LULC types. While, overall, the area of water was reduced by 59461.2 ha. The dominated changing direction of water body areas are vegetation and agricultural land. There were 124191.9 ha water body translated to water body and building area. The possible reason was because that reclaiming land and agricultural land from the river have reduced the area of the river that flows through the research area.(2)There was a growth-downward-growth trend of vegetation. During the period from 1996 to 2006, vegetation had a tendency of slight increase. And single dynamic degree of vegetation was 0.89%. Overall, the area of vegetation has decreased by 16,920.9 ha in the last 30 years. The main flow direction of vegetation was building area and bare land. And the conversion ratio was 4.31% and 2.85% for building area and bare land, respectively. The primary land types that translate to vegetation were building land, agricultural land, and bare land. The main reason were that urban growth and development are influenced by the national sustainable development strategy and the concept of green harmonious development. The reforestation activities make the area increasing. During the period from 2006 to 2016, the severe shortage of agricultural leading to vegetation transform into agricultural land.(3)Building land appeared a sustainable increasing trend. Building area of the research city increased by 1016336.7 ha from 1986 to 2016. The mainly transformation is that vegetation and farmland contributed more to building area than other LULC types. The primary reasons were that the urbanization process in the research city has been accelerating, and the urban expand was increasing by industrial infrastructure development. The expansion rate of urban building land is larger than reduced rate of building in villages, thus the urban building area increased overall.(4)There was a downward-growth-downward trend of farmland. Farmland mainly translated to vegetation and building land during this period. The conversion ratio of farmland was 14.57% from 1986 to 2016. And the main LULC types that translate to farmland were vegetation and bare land. That was because building land expansion occupied a large amount of farmland with the development of urbanization and industrialization. Furthermore, people in the village moved into the city, which cause some areas of farmland to be shelved and abandoned. The development of forestry reduced the area of farmland to a certain extent.(5)Bare land appeared downward-growth-downward trend. There was a sharp reduce of bare land from 1996 to 2006. The reduced rate of bare land was 3.54% during the period and was slower than the previous period. The total area of bare land increased by 56415.6 ha. The bare land was mainly converted into vegetation and farmland. The primary reasons were that the contradiction between population and resource is protruding, and the demand of urban land is increasing. These result in the exploit of unused land. Simultaneously, with the development of the economy and the afforestation from the bare land, the whole area of the bare land tends to decrease.

4.2. Prediction of Changing Trend

Based on the result of comparison between simulated land cover map in 2016 and the classification map in 2016, the prediction accuracy of the CA-Markov model was analyzed and verified (shown in Table 8). Results demonstrate that all relative errors between classified and predicted LULC type area are less than 8%. So it is scientific and reasonable to predict landscape changes based on the CA-Markov model. Then the CA-Markov model was selected to simulate land cover map in 2026 based on LULC maps captured in 2006 and 2016. The results are shown as Table 9 and Figure 3.

The results in Table 9 indicate that the largest changing amplitude of landscape type with a decrease trend is vegetation, which has an amount of 38.07 percent of total LULC types during the period of 2016 and 2026. The secondary change in the biggest of all is agriculture land, which will increase by 191006.1 ha from 2016 to 2026. The area of vegetation, farmland, bare land, and water body all will be increased at different degree based on the simulated results. It demonstrates that the urbanization in the research area is accelerating with the development of the urban economy. Simultaneously, a large reduction of vegetation shows that the ecological balance in the research area might be destroyed at a certain extent.

4.3. Landscape Pattern Analysis at Different Levels
4.3.1. Analysis of Landscape Pattern at Patch Level

Based on primary statistics, in total, the results shown an interesting phenomenon that relative large number of patches (NP) was bare land, the largest interspersion juxtaposition index (IJI), and the least edge density (ED) was water body. And the landscape types in the research, which has relative large largest patch index (LPI) was farmland and building area, while relative large aggregation index (AI) was vegetation and water body. The landscape types with relative large mean shape index (MSI) was vegetation, and the relative large fractal dimension index (FDI) was agricultural land and vegetation. So the landscape of various LULC types has changed significantly since 1986 shown from Figures 411 for NP, LPI, PD, ED, IJI, MSI, and AI, respectively. And the detailed analysis of landscape change can be discussed as follows: (1)The changing of farmland landscape. NP of farmland increased by 27,474, and patch density (PD) of farmland increased by 1.144 from 1986 to 2016. By prediction, largest patch index (LPI), aggregation index (AI), and interspersion juxtaposition index (IJI) will be reduced at different degrees by 2026. These results indicate that farmland is seriously disturbed by human activity in the process of urbanization and industrialization. Massive activities of economic construction and development occupied many farmlands. The increase of NP and the decrease of AI demonstrate that the patches were more fragmented, and the original patches were split up into discrete pieces of debris and distributed more widely. The transfer between farmland and other landscape types becomes more frequent. Dominance role of farmland is reducing, and landscape pattern of selected area was more and more broken. The degree of invading border is higher; land spatial structure will be more complicated.(2)The transformation of vegetation landscape. The mean shape index (MSI), fractal dimension index (FDI) and aggregation index (AI) of vegetation were relatively large from 1986 to 2016. At the meantime, AI and LPI of vegetation were decreasing during the period. Due to the ecological construction activities in the research area, the NP of vegetation increased and LPI was decreased. These show that human activities and interference were becoming more and more frequent in recent years. Human being efforts aggravated fragmentation of vegetation and increased the connectivity intensify. The fact that the rising of NP, IJI, and MSI of vegetation indicates that distribution tend of vegetation was disperse. The complex shape, means the deep irregular shape, which demonstrate human activity of developing vegetation land has been enhanced. According to results of prediction, NP and IJI will decline, AI will rise, and the vegetation landscape tends to be concentrated by 2026. The number of adjacent patches will increase, and the control power of the superior patches will decrease by that time.(3)The changing of bare land landscape. The bare land had the largest NP index and had the least LPI index during the period from 1986 to 2016. NP of bare land decreased by 24,402, ED of bare land reduced by 9.382, IJI declined, and AI increased. These indicated that the area of the bare land landscape was decreasing, and the distribution of bare land was more concentrated. According to the predicted results, NP of bare land will be further reduced and AI will rise significantly by 2026, which demonstrate that with the development of urbanization in the research area, the utilization of land will be become more efficienct and remarkable.(4)The landscape pattern changes of water body. Water body had the relative high AI and the largest IJI, while it has the least LPI and ED, during the period from 1986 to 2016. There was a sustained downward trend in NP and PD of water body, but there was a continued growth tendency in AI and IJI of water body. These demonstrates that due to interference of human activities, rainfall, urban construction, and other factors, the patch shape of the water body landscape tends to be simplified, the separated patches are joined together, and the distribution tend to be concentrated. The predicted results shown that NP will rise and AI and IJI will decline by 2026. These phenomena indicate that the area of the water body will increase, and the distribution become dispersed.(5)The changing of built-up landscape. In recent 30 years, the NP of building land increased by 39,155, and the FDI increase by 0.01. These numbers demonstrated that the area of building land is soar and occupying a domination because the increase in population and the development of urbanization. The phenomenon that AI and LPI of building land rise demonstrated the increasing connectivity for expansion of urban and the continuous distribution of large-scale area rural building land. The rising MSI of building land rising shown that the patches were clustered, the boundaries were fragmented, the distribution is more concentrated, and the shape is more complex. Based on the results of prediction, there will be a sustained downward trend in NP of building land, but a continued growth trend tendency in FDI of building land by 2026. These predicted results demonstrate that there is lack of reasonable plans for infrastructure construction in the preliminary stage, with the development of economic and growth of population, which makes the structure of patches more complex. While in the late stage, there is more scientific and reasonable construction plan, which makes the structure of patches not so much complicated.

4.3.2. Landscape Pattern Analysis at Landscape Level

In terms of the landscape view of the entire research area (shown in Table 10), there has similar changing trend of the total NP and FI, shown as downward-growth-downward trend. FDI has undergone a process from large to small and small to large, which is contrary to contagion index (CONTAG). The statistics shown that there is a sustained downward trend of LPI. The mentioned results above demonstrate that the dominate landscape in the previous stage has been replaced by the small patches gradually, and the landscape of the research area tend to be more complex and diversified. During the period from 1986 to 2006, the area of vegetation landscape increased significantly, the dominance of regional landscape was increasing, and the types of patches formed satisficed connectivity. While the construction of the research area lacks reasonable plans, the structure of the patches is more complex with economic development and population growth. But with the help of scientific and reasonable plan by decision maker, the structure of patches become simplification during the latter period (from 2006 to 2026). Simultaneously, with the rapid decrease of vegetation, the fragmentation degree of the landscape will increase, the landscape dominance will decrease, and the regional landscape will form a dense pattern with various elements.

During the first period from 1986 to 2006, there was a sustained downward trend in Shannon’s evenness index (SHEI) and Shannon’s diversity index (SHDI) shown that the heterogeneity of landscape decrease, the proportion of the dominant landscape rise, and the control function of the single component increase. The reason of the result is because that there is plain area domination in the selected research area, and that the policy of returning farmland to forest causes vegetation area increasing. The vegetation, as a dominant landscape type, whose proportion further improved and whose ability of landscape controlling has been enhanced significantly. During the second period (from 2006 to 2026), SHEI and SHDI rise observably. This indicate that interference intensity of human activities enhanced, the rare patches area increased, landscape diversity rises, the effect of controlling action of advantage components of landscape to overall landscape has decreased, the landscape stability decreased, and the degree of landscape heterogeneity is gradually increased.

The primary reasons were that to meet the objective requirements of urban expansion over the next decade, the economy has been developing rapidly over the past decade, and the greatly increasing building area occupied huge amount of vegetation and farmland that led to the increase of rare patches such as water body. On a whole, the proportion differences of different landscape were gradually narrowed, and the diversity of the landscape and the complexity of the structure increased significantly.

4.4. Analysis, Assessment, and Prediction of Ecological Security Status

With the help of the ecological security assessment system of above mentioned, maps of ecological security assessment classes of Zhengzhou city were obtained (shown as Figure 12).

Based on the primary results, the mean value of five phases of comprehensive security assessment from 1986 to 2026 were calculated (shown in Table 11), and the ecological security curve of Zhengzhou was drawn (Figure 13). The proportion of ecological security class of Zhengzhou city was statistical according to the assessment standard in Table 3; the changing situation is shown in Figure 14 and (shown in Table 12).

During the 40 years, the ecological security assessment value shows a “V” shape trend, decrease first and then increase. The ecological security assessment value of the research area was 7.01 in 1986, which means the ecological condition was at the “security” level. While the ecological security situation worsens significantly in 1996, and the ecological security was at “danger” level at the time. From 2006 to 2016, the ecological security evaluation value increased gradually, and the ecological security situation improved. According to the predicted results, the ecological security assessment value of the research area will further increase by 2026, but the ecological security situation will be still at the “criticality” level. The changing curve of ecological security shown that the ecological security value of the research area in 1986 is relatively gentle; the ecological security was more balanced. The two peaks of the ecological security value in the research area from 1996 to 2026 are approximately 5.5 and 8, which shown that the ecological security value of most areas is about 5.5 and 8. Based on the analysis of remotely sensed classification images and the maps of comprehensive assessment ecological security level, a conclusion can be drawn as follows. The area with low ecological security value concentrates mainly upon the water body and building land. The area whose ecological security value is 5.5 distributes in vegetation and farmland types, and the area whose ecological security value is around 8 concentrates in vegetation and bare land where human disturbance is weak.

The results of horizontal comparison of ecological security status shown that ecology of the research area is major in the “security” level in 1986, with a proportion of 38.81%. And the proportion of “morbidity” level of ecology is less than 1%. These conclusions indicated that a complete and stable ecosystem had established in the research area, but the ecological security areas with “morbidity” level are mainly focused in the northern of the research area due to the lack of reasonable management. In 1996, the “criticality” level ecology was dominant in the research area. And the percentage of “Morbidity” level ecology reached to 9.17%, exceeding the “security” level, which shown that the ecosystem stability of the research area is affected seriously by human disturbance. In 2006, the area with “criticality” level of ecology was more than 60%, the area with “security” and “stability” level of ecology was more than 10%. In 2016, the area with “criticality” level of ecology remained dominant in the research area. According to the predicted results, the area with “criticality” level of ecology will occupy a greater advantage by 2026, with an obviously trend of aggregating to “criticality” level.

Combining with maps of comprehensive assessment of ecological security level, a 40-year long-term vertical comparison analysis of ecological security assessment over the research area was carried out. Due to the rapid development of urbanization and the weak awareness of human ecological protection, during the period from 1986 to 1996, the “danger” level of ecology in the central and northern research area converted to “morbidity” level. At the meantime, the “stability” level of ecology transformed to “criticality” and “security” level. The ecological vitality of the research area decreased obviously, and the ecological security of the research area was deteriorated significantly. The rapid increase of population may be caused by rapid development of urbanization, the boom of economy, and the extensive construction of infrastructure. The increase of population pressure and the decrease of vegetation coverage led to the decline of ecological security in the research area overall.

During the period from 1996 to 2006, with the implementation of the sustainable development strategy, the increasing of vegetation coverage, the “morbidity” level of ecology changed to the “danger” and “criticality” level gradually. So, the ecological security situation has been alleviated in some extend. The primary reason is that with the promoting of strategy of “central rise” and sustainable development by Chinese administrate, the economic development model has changed from traditional economy to intensive economy, which strengthen effort of ecological and environmental protection significantly. Meanwhile, with the overall implementation of “one family, one child” policy, the contradictions between population and resources has been relieved in some extent. So, the landscape ecological security situation was improved.

From 2006 to 2016, the area with “danger” level of ecology increased slightly, and the area with “security” level of ecology reduced slightly, which indicated that with the further strengthening environmental protection by government, the improvement of awareness of environmental protection by residents in the process of urban development, ecological security has been further improved. It is predicted that by 2026, the “morbidity” level area will further reduce. In the meantime, the tendency of gathering to middle level of ecological security is more obvious with the increasing of “criticality” level. Therefore, it is an urgent problem to be solved that protecting the current security landscape, improving the landscape which was in “danger” and “morbidity” level.

5. Conclusions

In this work, the synthetically analysis of urban surface environment changes and landscape ecological indexes extends the application scope and applicable environment of the remotely sensed data. An ecological security evaluation system was established by collaborative processing of society, economy, nature, and environment factors. This study provides a first-time 40-year period, from past to future, comprehensive changing trend of typical LULC types over the selected area. The study proves that a growth-changing trend in building land and vegetation, a downward trend in water body, farmland, and bare land. During this LULC change process, the patch shape of the water body tends to be simplified, and the distribution of this landscape tends to be dispersed. The diversity of the landscape and the complexity of the structure had increased significantly overall. The ecological security evaluation value of Zhengzhou city in the last 40 years shows a “V” shape trend.

There are also some suggestions that can be drawn based on this research to realize the sustainable development of the research area. These suggestions include but not limit to the establishment of urban smart growth plan, strengthen of urban landscape ecological construction, controlling the total amount of urban land use scientifically, allocating the urban land area reasonably, following the conversion regulation of farmland to built-up land strictly, sustainably preventing and controlling soil erosion, optimizing the urban planning and design, and so on. The administrators are also recommended to promote the relationships between economic development and ecological balance by taking reasonable measures.

Conflicts of Interest

The authors declare no conflict of interest. The founding sponsors had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript; and in the decision to publish the results.

Authors’ Contributions

Pei Liu led the research work, proposed the idea and structure of this manuscript, wrote the instruction and conclusion parts, and integrated the contributions from all authors. Ruimei Han and Shoujun Jia integrated the modules and codes, performed the experiments, wrote LULC map classification and CA-Markov simulation parts and tested. Ruimei Han who implemented landscape indexes selection and extraction, established ESA system, wrote the corresponding parts, and also collaborated with Pei Liu to conceive and design the experiments. Hanwei Zhang was the coleader of the group and involved in the algorithm implementations through the paper, performed changing detection, and wrote the corresponding parts.

Acknowledgments

This work was supported by the National Natural Science Foundation of China (Grant nos. 41601450, 41601392, 41401403), the Key Research Project Plan of Colleges and Universities in Henan Province (Grant no. 16A420004), the Key Technology R&D Program of Henan Province (Grant nos. 182102310860, 162102310089), Henan postdoctoral foundation and Ph.D. fund of Henan Polytechnic University (Grant nos. B2015-20, B2014-18).