An improved urban cellular automata model by using the trend-adjusted neighborhood

Cellular automata (CA)-based models have been extensively used in urban sprawl modeling. Presently, most studies focused on the improvement of spatial representation in the modeling, with limited efforts for considering the temporal context of urban sprawl. In this paper, we developed a Logistic-Trend-CA model by proposing a trend-adjusted neighborhood as a weighting factor using the information of historical urban sprawl and integrating this factor in the commonly used Logistic-CA model. We applied the developed model in the Beijing-Tianjin-Hebei region of China and analyzed the model performance to the start year, the suitability surface, and the neighborhood size. Our results indicate the proposed Logistic-Trend-CA model outperforms the traditional Logistic-CA model significantly, resulting in about 18% and 14% improvements in modeling urban sprawl at medium (1 km) and fine (30 m) resolutions, respectively. The proposed Logistic-Trend-CA model is more suitable for urban sprawl modeling over a long temporal interval than the traditional Logistic-CA model. In addition, this new model is not sensitive to the suitability surface calibrated from different periods and spaces, and its performance decreases with the increase of the neighborhood size. The proposed model shows potential for modeling future urban sprawl spanning a long period at regional and global scales.


Introduction
Urban sprawl modeling is crucial for evaluating potential ecological and environmental risks caused by global urbanization in the future. By 2050, it is expected that the global urban population will reach~70% of the world's population (United Nations 2019). Such an explosive growth of urban population would result in the rapid expansion of urban extent, particularly in the fastdeveloping areas such as Asia and Africa. The rapid urban sprawl has many adverse effects on sustainable development, such as air pollution, agriculture land loss, deforestation, and public health (DeFries et al. 2010;Foley et al. 2005;Gong et al. 2012;Li et al. 2019b;Zhang et al. 2012). Therefore, to evaluate and minimize these risks and pursue cities' ecological and environmental sustainability, modeling future urban sprawl under different scenarios is highly required for policymakers and scientists to analyze the urban dynamics in complex resource-constrained environments and then make good decisions for city planning and management Li and Gong 2016).
Cellular automata (CA)-based models have been widely explored in urban sprawl modeling for simplicity, transparency, and flexibility (Santé et al. 2010). The key of the CA-based urban sprawl model (i.e., urban CA model) is the self-evolution of the urban cell that is driven by its neighbors (Batty and Xie 1994). Thus far, a number of CA-based urban sprawl models have been developed and continuously improved through extensive studies, such as the SLEUTH (Slope, Land use, Exclusion, Urban, Transportation, and Hillshade) model (Clarke et al. 1997), the Logistic-CA model (Wu 2002), the Constrained CA (Li and Yeh 2000), the Fuzzy CA (Liu and Phinn 2003), and the Markov CA (Shafizadeh Moghadam and Helbich 2013). By virtue of different approaches to represent the driving factor of urban development by transition rules and neighborhoods, these CA-based urban sprawl models can produce an urban morphology that is close to actual urban lands (Li et al. 2017b;Liao et al. 2019). Among these models, the Logistic-CA model is popular because of its easy implementation and explanation.
Transition rules and neighborhoods are two critical components in urban CA models. The original transition rules of urban CA models are explicit "if-then" regulations. These rules were later developed as probabilitybased rules (or called suitability surface) that integrate multiple spatial proxies (Li and Gong 2016;Santé et al. 2010). Different approaches and models, including empirical methods such as regressions and multi-criterion estimation (Hu and Lo 2007;Wu 1998Wu , 2002 and nonlinear models such as decision tree or artificial neural networks (Li and Yeh 2004), have been widely explored to obtain the suitability surface. Although non-linear models (e.g., machine learning approaches) can achieve a better performance than the empirical regression approach, the resulting suitability surface is difficult to explain. For the neighborhood of urban CA models, most studies focus on its spatial configuration (Kocabas and Dragicevic 2006). The shape and size are two common indicators used to form the configuration of the neighborhood. Additional information within the neighborhood, such as the distance to the central cell (van Vliet et al. 2009) and the land use/cover composition (Wu et al. 2012), was also widely studied in the urban CA modeling.
At present, limited efforts have been made to explore the temporal trend of urban sprawl, although many urban CA models have been developed. These developed urban CA models mostly focused on improving the model capacity regarding the spatial information, e.g., the distance of the central cell to surrounding urban infrastructures (Santé et al. 2010;van Vliet et al. 2009). Despite the growth of cities generally follows its historical pathway, the temporal context of urban development was seldom included in the urban CA models. That is, compared with regions developed at early years, the more recently developed regions have a higher likelihood to be developed in the near future (Liu et al. 2017). Satellite observations also confirmed this trend from a long temporal perspective at national and global scales (Gong et al. 2019;Huang et al. 2020;Li et al. 2015a;Zhou et al. 2018;Zhou et al. 2014;Zhou et al. 2015). Thus, the information of the temporal context of urban development shows a great potential to improve the performance of urban CA models.
Motivated by this idea, in this paper, we firstly developed an improved urban CA model using a trendadjusted neighborhood, of which the historical pathway of urban sprawl was considered. We applied this newly developed urban CA model in the Beijing-Tianjin-Hebei region, a rapidly developing region in the North Plain of China.

Study area
The Beijing-Tianjin-Hebei region, in the north of North China Plain, is the largest urban metropolitan areas of China ( Fig. 1). This region occupies an administrative area of 216,600 km 2 (Dong et al. 2008), with more than 100 million people. During the past decades, this region experienced rapid growth of population and economy, resulting in a notable sprawl of urban extent. As two primary engines of this region, Beijing and Tianjin lead the development in this area with almost an exponential growth of urban areas over past decades (Chai and Li 2018;Li et al. 2015a). This rapid urbanization is raising public concerns on water scarcity (Li et al. 2018a), energy consumption (Wang and Chen 2016), and air pollution (Liu et al. 2018). Therefore, modeling of historical urban sprawl and predicting of future growth are urgently needed in this region.

Data collection
We collected seven spatial proxies in the Beijing-Tianjin-Hebei region in the urban sprawl modeling (Table S1). These proxies consist of spatial features and images such as terrain, land cover, traffic, and location. Except for proxies in the category of terrain (i.e., elevation and slope), all other proxies were processed as the distance to city centers, roads, and specific land cover types (Fig. 2). These spatial proxies were used to determine the suitability of urban sprawl of each pixel, according to their biophysical and socioeconomic conditions .
We used the annual urban extent data derived from nighttime light (NTL) and Landsat satellite data in the model calibration and evaluation. These two annual urban extent datasets were used not only for providing the temporal trend of urban sprawl but also for evaluating the robustness of the improved urban CA model. The NTL-derived urban extent maps (Fig. 2h) span from 1992 to 2013, with a medium resolution of 1 km. The mean accuracy of this developed dataset is about 89% in China (Zhou et al. 2018). The Landsat-derived urban  , distance to city centers (c), distance to highways (d), distance to major roads (e), distance to local roads (f), land cover (g), and urban extent dynamics from NTL observations (h)

Li et al. Ecological Processes
(2020) 9:28 extent maps (Fig. S1) have a longer temporal interval  than the NTL-derived results, with a fine spatial resolution of 30 m. Based on the long-term Landsat time series data, this Landsat-derived dataset was generated by a temporal segmentation approach (Li et al. 2018b). The overall accuracy of the detected urbanized year in this dataset is around 83% in the Beijing-Tianjin-Hebei region. Besides, both of these two urban dynamic datasets follow the logic of urban development.
That is, this development is a monotonic conversion from non-urban to urban (Li et al. 2015a).

Framework
We developed a Logistic-Trend-CA model and assessed its performance to relevant factors in urban CA model ( Fig. 3). First, we proposed a trend-adjusted neighborhood with the consideration of the historical pathway of urban sprawl and developed a Logistic-Trend-CA model. Second, we analyzed the performance of the proposed model to key elements in the urban CA model, including the start year of modeling, the suitability surface, and the neighborhood size. The improved urban CA model was applied at 1-km and 30-m spatial resolutions to explore its capability in cross-scale modeling. Details of each step are given in the following sections.

The logistic-trend-CA model
The urban CA model is a grid-based self-evolution system to simulate the dynamics of urban land (Batty and Xie 1999). In this system, the status (i.e., urban and nonurban) of each grid is determined by its surrounding neighbors. A non-urban grid is more likely to change to urban in the near future if there are more urban grids surrounded. Evolution of massive grids using this rule simultaneously can simulate the change of complex urban landscapes. With the consideration of additional spatial factors such as traffic networks and land covers, the urban CA model can be used to simulate the dynamic of urban land with a high degree of reliability. We built our Logistic-Trend-CA model on the Logistic-CA model as it has been widely used in urban sprawl modeling due to its clear explanation of spatial proxies and ease of implementation (Hu and Lo 2007;Wu 2002). The logistic regression function is the key of the Logistic-CA model. Its output is a spatially explicit suitability surface, which indicates the suitability for development under considerations of different spatial proxies. Assuming there are n spatial proxies [x 1 , x 2 , … x n ], the logistic regression function can be expressed as Eqs.
where P suit is the obtained suitability of development from the biophysical and socioeconomic conditions and b i and x i are the ith coefficient and spatial proxy, respectively. We improved the neighborhood of the Logistic-CA model by considering the historical pathway of urban sprawl. The neighborhood is a crucial component in the urban CA model because it is a basic driver of urban dynamics modeling (Kocabas and Dragicevic 2006). The configuration of the neighborhood is closely related to its size, shape, and surrounding land cover types. Here, we developed the trend-adjusted neighborhood by incorporating the historical pathway of urban sprawl as a weighting factor, based on the widely used Moore configuration (Eqs. (3 and 4)).
where Ω is the influence of neighborhood that considers the historical pathway of urban sprawl using a weighting factor of W ts ij . N u ij is the accumulated year of cell (i, j) with the status as urban from the annual urban time series data with a temporal interval of N. As a result, for a potential cell, urban neighbors that were developed in more recent years have larger impacts than those developed in earlier years. m is the window size, and Con() is a conditional function, which returns 1 when the status of cell (i, j) is urban.
Compared with traditional neighborhoods, the developed trend-adjusted neighborhood can result in a sprawl pattern following the historical pathway, as illustrated in Fig. 4. Urban sprawl has an inertia of development as it generally follows the temporal trend of historical development (Liu et al. 2017), i.e., there is a relatively higher development probability around those newly developed urban areas. As a result, the weighting factors of W ts ij of urban pixels developed in more recent years are higher than those developed in earlier years. Assuming the weighting factors of urbanized pixels across years are different as illustrated in Fig. 4 (a), thus, pixels 1 and 2 have the same neighborhood influence if using the traditional neighborhood. However, if taking into account of the historical pathway of urban sprawl, pixel 2 has a higher neighborhood influence because its surrounding pixels were developed more recent compared to pixel 1 ( Fig. 4 (b)), although they have the same number of urban neighbors for pixels 1 and 2. Such smaller difference regarding the neighborhood influence would result in a notable different sprawl pattern after several iterations (Fig. 4 (c)), as most urbanized pixels were developed following the historical pathway if the neighborhood was adjusted by the temporal trend.
We also included land constraint and stochastic perturbation in the developed Logistic-Trend-CA model. Restricted lands, such as water and protected areas, were not considered for development in our model; thus, they were represented as a land constraint term as Land = 0 ). In addition, we used the stochastic perturbation SP to represent unconsidered factors (e.g., policy) in the modeling (White and Engelen 1993), as expressed in Eq. (5).
where SP is the stochastic perturbation, λ is a random value [0, 1], and α is a parameter to determine the degree of perturbation.
The development probability was calculated based on the suitability surface, neighborhood, land constrain, and stochastic perturbation. For urban time series data derived from NTL and Landsat, we determined their Fig. 4 Illustration of the trend-adjusted neighborhood in the urban CA model. Urban pixels with different developing years and weights (a), calculation of neighborhood impacts for pixels 1 and 2 (b), and modeled urban sprawl using different neighborhoods (c). This illustrative figure shows the results of allocating three urban pixels per iteration, with a total of nine urban pixels. Trend weights in (a) were calculated using Eq. (3) with a temporal span of 4 years development probabilities using Eqs. (6) and (7), respectively. The SP is not considered in the modeling with NTLderived urban extent maps, due to the homogeneity of urban land within the boundary (Zhou et al. 2018).

Model validation
We assessed the performance of the developed Logistic-Trend-CA model to key factors in urban CA model. Two quantitative metrics were used for the assessment, namely the overall accuracy (OA) and the figure of merit (FOM). The OA was directly calculated as the percentage of consistent pixels to all pixels in the entire map, while the FOM indicates the consistency between modeled and observed maps on those changed pixels. The FOM has been widely used in many studies of urban CA models since it can provide a relatively comprehensive evaluation of the model performance (Chen et al. 2014;Li et al. 2014;Pontius et al. 2007). The FOM can be expressed as Eq. (8) (Pontius et al. 2008).
where FOM is the figure of merit, B is the number of observed urban pixels that were simulated as urban, A is the number of observed urban pixels that were simulated as non-urban, and C is the number of observed non-urban pixels that were simulated as urban.
We evaluated the model performance by exploring sensitivities of three key factors in our urban CA model: the start year of modeling, the suitability surface, and the neighborhood size. The influences of these three key factors can be quantitatively evaluated with clear physical meanings in the CA model. Although there are other factors that may also influence the model performance, they are considered in the selected elements in the model. For example, the urban spatial configuration can be captured by the neighborhood and stochastic disturbance in the CA model. Other factors (e.g., the degree of urban development) are more related to regional urban demand compared to the spatial allocation of increased urban demand (Li et al. 2019a). First, we examined the modeling capability of the Logistic-Trend-CA model over a long temporal span through changing the start year of modeling, which is closely related to the iterations and error propagation in the modeling ). Second, we investigated the model performance to suitability surfaces using calibrated results from different periods. The suitability surface characterizes the likelihood of urban development from the biogeophysical (e.g., terrain and land cover) and socioeconomic (e.g., traffic networks) aspects. In addition, we evaluated the derived suitability surface using the receiver operating characteristic (ROC) approach (Liu et al. 2017;Pontius Jr et al. 2001;Wu et al. 2009). The ROC curve was calculated by dividing the continuous suitability surface into binary maps using different thresholds and then comparing the derived binary map with the reference map. Finally, we compared the model performance by varying neighborhood sizes. The neighborhood is a crucial component to drive the self-evolution of urban land system, and the neighborhood size reflects the degree of local impact from neighbors (Kocabas and Dragicevic 2006). Many urban CA models have been developed for particular applications with different structures, functions, and data requirements (Li and Gong 2016). Quantitative indicators such as the FOM have been used for comparing urban CA models. We evaluated our model performance based on the FOM and compared it with previous studies. In addition, we compared our model with the similar Logistic-CA model, which has been widely used in previous urban CA studies and can serve as a benchmarking model, for several key factors. The Logistic-CA model also has the same structure as our proposed model except for the consideration of the neighborhood.

Setting of the urban CA model
The inputs of urban CA model are the urban extent map in the beginning year associated with a variety of spatial proxies (Fig. 2) and a set of parameters (Li et al. 2017a), and the output of our model is the urban extent map in the target year. In our study, the neighborhood size (m) was set as 3 and 5 in calculating the influence of neighborhood (Ω) using Eq. (4), for urban sprawl modeling with medium (1 km) and fine (30 m) resolutions, respectively. The degree of stochastic perturbation (α) was set as 3 as suggested for modeling at a 30-m resolution using Eq. (5) . Also, we set the restricted conversion type as water in our study to avoid the conversion from water to urban. Finally, the modeled results were compared with the observed urban extent map from remote sensing observations in the same year.

Results and discussion
Performance of the logistic-trend-CA model The developed Logistic-Trend-CA model using the trend-adjusted neighborhood outperformed the traditional Logistic-CA model. Improvements of FOM are about 18% and 14% for urban sprawl modeling at medium (1 km) and fine (30 m) resolutions, respectively (Figs. 5 and 6). The OA increases by around 2-3% using the Logistic-Trend-CA model. This suggests the developed model is robust regarding the model performance at different resolutions. Specifically, the omission (i.e., pixels observed as urban but simulated as non-urban) and commission (i.e., pixels observed as non-urban but simulated as urban) errors were considerably reduced when the historical pathway of urban sprawl was considered (Figs. 5 and 6). That is, the trend-adjusted neighborhood improves the urban sprawl pathway during the modeling, which further reduces the error generation and propagation Yeh and Li 2006), particularly for modeling over a relatively long temporal interval. As a result, the developed Logistic-Trend-CA model can simulate more realistic urban forms (or landscapes) compared to the traditional Logistic-Trend-CA model. However, it is worth to note that the improvement of Logistic-Trend-CA model is related to the time span of the used temporal information in the neighborhood. The improvement of FOM in Fig. 5 would decrease to 10% and 6% when the time span reduces to 10 and 5 years, respectively. If using a very limited time span (e.g., 1-2 years), outputs from Logistic-Trend-CA and Logistic-CA models are almost the same.

Validation of model performance to key factors
The start year of modeling A shorter modeling period with the start year of 2002 can increase the OA but decrease the FOM of the modeled result compared to a longer modeling period with the start year of 1992. Since the time step of our modeling is annual, the impact of start year on the model performance is mainly determined by the number of iterations, which further affect the error generation and propagation in the modeling . The comparison of model performance using different temporal intervals indicates an opposite trend of OA and FOM (Fig. 7 and Fig. S2). That is, the OA is slightly increased from 97 to 98% while the FOM is decreased from 52 to 44% when using a relatively shorter modeling period at a medium resolution (1 km). Similarly, for urban extent maps with a fine resolution (30 m), the improvement of OA is about 2% while the decrease of FOM is about 8-9%. Such a contrastive trend of OA and FOM is related to the consistency of urban extent maps between the start and end years. For example, although the FOM of modeling with a start year of 2002 ( Fig. 7b) is lower than that with a start year of 1992 (Fig. 7a), the initial urban extent in 2002 has excluded errors generated and propagated from 1992 to 2002, resulting in a higher overall agreement.
The improvement of the Logistic-Trend-CA model relative to the Logistic-CA model decreases when modeling with a relatively short period (Fig. 8). Introducing the historical pathway of urban sprawl improves the performance using the Logistic-Trend-CA  (Fig. 8a). Similar conclusions are also confirmed in modeling cases of Beijing and Tianjin at the 30-m spatial resolution during the modeling period of 2005-2015, i.e., the FOM is increased from 53 to 54% in Beijing and from 32 to 39% in Tianjin (Fig. 8b, c). Although the improvement of FOM is reduced by about 8% when the modeling period is shortened from 1992-2013 (Fig.  5) to 2002-2013 (Fig. 8a) in the Beijing-Tianjin-Hebei region, the Logistic-Trend-CA model remains a good performance compared to the Logistic-CA model, particularly for modeling cases with relatively worse performances of the suitability surface (i.e., Tianjin) (Fig. 9c). Also, it should be noted that the performance of the proposed model is related to urban expansion patterns in cities, and the improvement of the Logistic-Trend-CA model for a concentric growth across different directions around the urban center is not significant (e.g., Fig. 8b).

The suitability surface
Our Logistic-Trend-CA model is not sensitive to suitability surfaces that are derived during different periods. These suitability surfaces were calibrated using the logistic regression model, based on training samples collected from urbanized and persistent regions in different periods. We found the temporal effect of training samples collected from different periods on the derived suitability surface is limited (Fig. 9), i.e., their ROC curves are close. The impact of different suitability surfaces on urban sprawl modeling is  (Table 1), i.e., the FOMs are similar in these experiments. Modeled results at both the 1-km and 30-m spatial resolutions are almost the same using suitability surfaces in different periods. This is because the development probability of urban sprawl in each iteration is jointly determined by the neighborhood and the suitability surface. The suitability surface is assumed to be persistent, and its contribution is reduced as errors propagated during the modeling (Li et al. 2015b;Santé et al. 2010), whereas the trend-adjusted neighborhood is updated iteratively and plays a dominant role in the modeling. It should be noted that suitability surfaces used in our study were mainly derived from the logistic regression model, which is a statistical model and is more robust for training samples compared with data miningbased approaches such as random forest and neural networks ).

The neighborhood size
The performance of Logistic-Trend-CA decreases with the increase of the neighborhood size (Fig. 10). This finding is consistent for urban sprawl modeling at different spatial resolutions (i.e., 1 km and 30 m). The FOM is highest when the neighborhood size is 3; thereafter, it decreases with the increase of the neighborhood size. This relationship is more distinctive when the window size is small (e.g., lower than 7), suggesting the contribution of included weighting factors from the historical pathway of urban sprawl decreased when increasing the neighborhood size (Eq. 4). Accordingly, the improvement of the Logistic-Trend-CA model relative to the Logistic-CA model was reduced. In addition, urban sprawl modeling is more sensitive to the neighborhood size with a relatively coarse spatial resolution, comparing modeled results at the 1-km (Fig. 10a) and 30-m (Fig. 10b, c) spatial resolutions. The impact of local neighbors on the model performance decreases due to the decrease of neighborhood intensity of the central pixel, with the increase of neighborhood size.

Conclusions
In this study, we developed a Logistic-Trend-CA model with the consideration of the historical pathway of urban sprawl and tested it in the Beijing-Tianjin-Hebei region of China. In this model, we proposed a trend-adjusted neighborhood as a weighting factor using the historical pathway of urban sprawl. This improved neighborhood was integrated with the widely used logistic regression function to simulate urban sprawl. We applied this model in the Beijing-Tianjin-Hebei region using the time   Performance (i.e., ROC curves) of derived suitability surfaces using training samples collected in different periods in the Beijing-Tianjin-Hebei region at the 1-km spatial resolution (a) and in Beijing (b) and Tianjin (c) at the 30-m spatial resolution. The suitability surface with a ROC curve in the upper-left corner indicates a good performance series data of urban extent from NTL and Landsat observations. Model performance was evaluated and compared with the traditional Logistic-CA model. The robustness was explored through analyzing the model performance to key factors in urban CA model.
We found our Logistic-Trend-CA model notably outperforms the traditional Logistic-CA model. The improvement of FOM is around 18% and 14% using the Logistic-Trend-CA model at the 1-km and 30-m spatial resolutions, respectively, compared to the traditional Logistic-CA model. The Logistic-Trend-CA model performs well for modeling studies with a long temporal span. In addition, it is not sensitive to suitability surfaces derived from the logistic regression model in different periods, and the trend-adjusted neighborhood plays an important role in the modeling. Finally, the performance of the Logistic-Trend-CA model decreases with the increase of neighborhood size.
This study opens a new research avenue to incorporate the temporal context information in urban CA models. Through using the temporal context information of historical urban sprawl, the Logistic-Trend-CA model shows a good performance in simulating future urban sprawl with a long interval (e.g., decades). Also, the developed urban CA model in this study performs well at the 1-km spatial resolution, showing its capabilities in global urban sprawl modeling using the time series data of urban extent from NTL observations (Zhou et al. 2018). Thus, our developed urban CA model can be used for urban expansion modeling with a long temporal span because historical growth of urban extent can be incorporated into the modeling. Such improvement can mitigate the uncertainty in modeling urban growth using the information of temporal contexts. However, it is worth to note that the temporal effect of suitability surface (e.g., road expansion and land cover change during the modeling period) is not considered in our model, although we included the dynamics of urban extent in the neighborhood component. Also, our model needs improvement for simulating policy-induced changes in urbanized areas without neighboring urban pixels (Liu et al. 2010), which could occur in rapidly developing regions in China. The corresponding improvements for these common limitations in CA models are needed for our Logistic-Trend-CA model in future studies.