Patch-Level and Neighborhood-Dependency Spatial Optimization Method (PNO): Application to Urban Land-Use Planning to Facilitate Both Socio-Economic and Environmental Development in Beijing

Abstract

Rapid urban expansion and chaotic urban land-use patterns cause many socio-economic and environmental issues, e.g., traffic congestion and urban heat islands; thus, scientific planning considering land-use trade-offs and layout optimization is highly required for resolving these issues, especially in the urban renewal stage. However, previous spatial optimization methods were weak in processing land-use patches and ignored their neighborhood dependency, leading to fragmented and inapplicable optimization results. Accordingly, this study proposes a patch-level and neighborhood-dependency spatial optimization method (PNO) to adjust urban land-use patterns considering multiple optimization targets (i.e., improving population and economy but controlling land surface temperature). The PNO represents land-use patterns in a graph structure, quantifies land-use patterns’ impacts on the population, economy, and land surface temperature, defines the spatiotemporal constraints of land-use optimization considering neighborhood-dependency and optimization sequences, and finally optimizes land uses and their spatial layouts based on a multi-objective genetic algorithm. Experiments were conducted in the urban area of Beijing, and the results suggested that, after optimization, the population and GDP can be improved by 667,323 people (4.72%) and USD 10.69 billion in products (2.75%) in the study area; meanwhile, the land surface temperature can be reduced by 0.12 °C (−0.32%). Through comparison, the proposed PNO outperforms previous spatial optimization methods, e.g., NSGA-II, in processing land-use patches as well as their neighborhoods. Taking the land-use map in 2022 as a reference, the PNO optimization results are more consistent with actual land-use changes (consistency of 25%), compared to the existing spatial optimization results (consistency of 10.6%). Thus, PNO is more applicable to land-use planning in urban renewal circumstances.

1. Introduction

After decades of rapid urban expansion, many countries, particularly in the Global North, are now undergoing urban renewal, which fundamentally aims to trade off and optimize land-use patterns (LUPs) for both improving urban socio-economic carrying capacities and protecting local environments, playing a critical role in achieving sustainable urban development [1].

LUPs refer to land-use categories and their spatial layouts. During periods of rapid urban expansion, LUPs often become chaotic and disordered, as different land-use patches are usually developed in different stages and have different targets either for socio-economic development or ecological environmental protection. However, unscientific LUP development can result in various issues. For example, a high proportion of commercial land may stimulate economic growth but can also result in significant traffic congestion [2], and as cities develop and expand, jobs–housing spatial mismatches tend to worsen, leading to longer commutes, increased traffic congestion, and higher levels of air pollution [3]. Urban structures trap heat and accumulate greenhouse gases, aggravating the urban heat island effect and increasing the land surface temperature (LST) due to a lack of green spaces [4,5]. Urban renewal aims at adjusting the previous LUPs to address these issues to promote high-quality and sustainable development. Previous methods merely converted the land-use attributes pixel by pixel or patch by patch. However, these methods did not regard the LUPs as an integrated system and ignored the inner connections between the patches. Differently, this study proposes a patch-level and neighborhood-dependency spatial optimization method, which conducts overall adjustments to LUPs for a group of patches, including the land-use attributes of patches and relationship between patches (Figure 1).

How can we consider multiple targets (e.g., socio-economic development and environmental protection) and the current spatial layouts to trade off diverse land uses and optimize LUPs for urban renewal? It is a wide concern but also a difficult NP-hard issue [6,7]; thus, LUP optimization has become a hot topic all over the world. For example, early research about LUP optimization is a linear programming model first introduced by Schlager [8]. This model is simply target-oriented and not applicable to the nonlinear land-use optimization problem because the objectives of LUP optimization are competing and even, sometimes, conflicting [9]. To address the conflicts between the objectives in multi-objective optimization problems, some heuristic optimization algorithms were proposed and applied to the LUP optimization problem, such as the Ant Colony Algorithm (ACA), Particle Swarm Optimization (PSO), and the genetic algorithm (GA), which allow a trade-off between different objectives and can be applied to find near-optimal solutions [10]. The most widely used is the GA, which simulates the crossover operation and mutation operation in natural selection, evaluates the advantages and disadvantages of individuals by using the fitness function, and selects the individuals with higher fitness according to the principle of the survival of the fittest. In decision-making processes, the weighted-sum approach and Pareto-based methods are two commonly used techniques for multi-objective optimization. The weighted-sum approach is simple and computationally efficient, as it converts multiple objectives into a single target value by the weights given to each normalization objective value. Cao et al. [11] came up with a multi-objective land-use optimization model using a boundary-based genetic algorithm and employed the weighted-sum approach to select excellent land-use scenarios. However, this approach is highly sensitive to and dependent on the selection of weights and could lead to highly biased solutions [12]. Compared with the weighted-sum approach, the Pareto-based approach can realize a reasonable trade-off and balance among different optimization objectives and effectively avoid the biased solution. On the basis of this, Deb et al. [13] proposed NSGA-II, where individuals are classified and sorted by several Pareto sets through rapid non-dominated sort, and any individual in the same Pareto set cannot dominate another individual completely, which means that any individual cannot be more optimal than any other individual in all aspects [14]. By using rapid non-dominated sort, the time complexity of the genetic algorithm is from O(n³) to O(n²). NSGA-II also provides an effective solution for nonlinear parameter regression and parameter optimization problems [15]. Based on NSGA-II, many complicated and conflicting multi-objective issues in various fields can be solved easily: Wang [16] proposed a framework of support vector regression NSGA-II to reduce energy consumption in buildings on the premise of ensuring indoor thermal comfort; Thakur [17] developed a multi-path refrigerated fresh multi-item eco-routing model using NSGA-II with hybrid chromosomes to maximize profit and minimize carbon emissions; Liu [18] put forward a low-carbon-oriented land-use allocation optimization model based on NSGA-II and generated a series of optimal patterns to support decision making for the low-carbon adjustment of the land-use spatial pattern; and Luan [19] took both ecological and economic benefits into account and then optimized land use for the Liangjiang New Area, China, based on NSGA-II.

As demonstrated above, NSGA-II can effectively balance multiple conflicting objectives and proves to be a powerful tool for achieving comprehensive and balanced optimization in scenarios. However, from the perspective of the structure and the inner relationship of LUPs, existing optimization methods are weak in shifting LUPs, as they cannot process land-use patches considering their neighborhoods and conversion rules. Accordingly, a patch-level and neighborhood-dependency optimization method is highly required but has three technical issues:

• (1). The urban land-use system is an integrated system, so the pattern may change a lot even if a patch changes. For example, assuming that a patch was changed to a transportation hub, its surrounding pattern would change correspondingly, and there would be more traffic in the neighborhood. Thus, a proper optimization method contraposing LUP optimization based on patches should be developed.

• (2). The relationship between optimization objectives and LUPs is unclear and hard to quantify; thus, LUP optimization results’ impacts on the socio-economic and environmental features are difficult to measure and simulate. Under this circumstance, a quantification method should be developed for representing LUPs and their relationship with multiple optimization objectives, e.g., population, gross domestic product (GDP), and LST in this study.

• (3). Current spatial optimization methods are usually designed based on pixels [20] or on a patch-by-patch basis without considering spatial dependencies between the neighboring patches [21], which does not represent actual urban land-use planning or renewal. In addition, there can be a lot of constraints for LUP optimization, for example, it is unacceptable to take any patch into a shanty town, which however have not been considered by current spatial optimization methods. Accordingly, a novel spatial optimization method needs development at patch level and considering diverse spatiotemporal constraints.

In summary, LUP optimization is highly difficult, and previous ACA, PSO, and GA optimization methods cannot resolve it. Accordingly, this work proposes a patch-level and neighborhood-dependency spatial optimization method (PNO) which generally has two characteristics: (1) PNO processes land-use patches considering their spatial patterns and thus can improve the population and GDP in the limited urban land, also reducing the LST; (2) the PNO results outperform previous optimization methods, e.g., NSGA-II, with more complete and applicable optimization results for supporting urban renewal; and (3) with reference to the land-use map in 2022 and the land-use changes from 2020 to 2022, the PNO results are more consistent with the actual land-use changes (25% consistency), which are more than double that of the existing spatial optimization results (10.6% consistency). To the best of our knowledge, this is the first method that measures neighboring relationships for land-use patch optimization, and it significantly contributes to urban land-use planning and renewal.

2. Materials and Methods

2.1. Research Area

Beijing, the capital of the People’s Republic of China, is the political, economic, and cultural heart of the nation. With a population approaching 21.89 million, Beijing has a land area of approximately 16,400 km² (Figure 2). Dominated by a temperate monsoon climate with four distinct seasons, Beijing has hot summers. In 2023, the temperature surpassed 35 °C for a total of 30 days. Over the past 60 years of rapid urbanization, the mean daily air temperature has risen by 4.4 °C, and the city is confronted by urban heat island effects, which are expected to be more and more severe over the coming years.

2.2. Data Sources

This study conducts an investigation within the 5th ring road in Beijing which covers about 735 km² and carries over 10 million persons (Figure 2). To optimize the land use across the study area, three geographic datasets including point-of-interest (POI), land-use patch, and road data are taken into consideration.

A satellite image (Figure 2) was obtained from the Gaofen-1 satellite in 2020, covering the entire study area with a high resolution of 2 m and including three visible bands (i.e., red, green, and blue). Furthermore, the image had been ortho-rectified.

Land-use data from 2020 (Figure 3a) and 2022 (Figure 3b) were used. The land-use data consisted of different categories of polygons, including forests, water, unused lands, transportation, parks, industrial lands, institutions, commercial lands, farmland, and first-, second-, and third-level residential districts. The data are the result of a multimodal semantic segmentation proposed by Zhang et al. [22]. The overall classification accuracy of this model is 90.8%

Roads (Figure 3c) consisted of polylines collected from AutoNavi Map in May 2018. The roads are related to urban functions closely, and road network density is associated with population density, GDP, and other socio-economic factors. What is more, artificial surfaces including roads have poor water permeability and low specific heat capacity and thereby can cause high surface temperatures. Thus, roads and population density, GDP, and LST are tightly relevant.

The population data (Figure 3d) are a raster image covering the entire study area which was acquired from the website worldpop.org in 2020. The pixel value represents the total number of persons in each pixel.

The GDP data (Figure 3e) are a raster image covering the whole study area which is the result of current research conducted by Zhao et al. [23]. The pixel value represents the gross domestic product (billion dollars) in each pixel which was forecasted by using nighttime light time series and population images.

The LST data (Figure 3f) are remote sensing synthetic data and one product of the global land surface satellite (GLASS) dataset in 2020, which is a high-resolution, high-precision, and long-term (over 35 years) global surface remote sensing product. The product is generated through inversion algorithms using multi-source remote sensing data and ground measurement data. The pixel value represents the LST in each pixel.

The POI data (Figure 4) contain 253,253 POI points collected from Baidu Map in 2020, providing timely updates to POI data. Thus, they can guarantee that the functional categories and locations of POIs are precise and dependable. As shown in Figure 4, these POIs are reclassified into 4 categories, consisting of Commercial House, Enterprises, Medical Service, and Science, Education, and Culture.

3. Research Methods

This study essentially aims to propose a new method based on patch data, which uses a Non-dominated Sorting Genetic Algorithm (NSGA-II) and regards population, GDP, and LST in summer as optimization objectives to realize the multi-objective optimization. Meanwhile, this method takes the impact of crossover and mutation sites on their neighbor patches into consideration during the process of crossover and mutation. This section mainly contains 5 parts, including the representation method of LUP, multi-objective simulation based on LUPs, the patch-level and neighborhood-dependency spatial optimization method, the spatiotemporal constraints of LUP optimization, and the effectiveness evaluation of PNO results (Figure 5).

3.1. Representation Method of LUPs

3.1.1. Category System of Land Uses

The study area is represented by several different polygons, each of which corresponds to a record which contains some attributes, such as the ID, area, and land-use category of this polygon. In this study, we only transform the land-use category of the patches, with their location, shape, size, and neighborhood relationships remaining unchanged. In order to represent LUPs in the GA, a chromosome is assumed as a pattern of land-use planning. The ID of each gene corresponds to the ID of a land patch, and the value of each gene can represent the land-use category of the patch (Figure 5a). Thus, the mathematical form of the land-use pattern in the study area can be expressed as follows:

(1)$Land=\left[\begin{array}{cccc}l{u}_1& l{u}_2& \dots & l{u}_n\end{array}\right]$

where n is the ID of each feature. lu_n is the code representing the land-use category of the feature whose ID is n. According to the “code for classification of urban land use and planning standards of development land” issued by the Ministry of Housing and Urban-Rural Development of the People’s Republic of China, this study considers 12 categories of land uses. Table 1 lists the land-use categories for each code and their descriptions. Among these categories, only 6 categories (“Parks”, “Industrial lands”, “Commercial lands”, “L1-res”, “L2-res”, “L3-res”) of patches are selected for optimization, and we only change the attribute of the patches from one type to one of the other five types in this study.

3.1.2. Graph Structure for Presenting LUPs

Graphs are a universal data format and an effective tool for modeling and analyzing intricate relationships among entities or objects consisting of vertices and edges, where the vertex characteristics represent different objects, and the edges that connect two nodes can describe complex relationships between objects. This study uses the undirected graph structure to present the LUPs and the neighborhood relationships of patches (Figure 5a). And each vertex represents the ID of the corresponding patch, while the edges between two vertices represent the neighborhood relationship of two patches. If two vertices are connected by an edge, it indicates that the two patches are adjacent. Taking a local area in the study area as an example, Figure 6 presents how the adjacencies of patches in LUPs (including A, B, C, D, E, F, G, H, I, and J) are represented using an undirected graph structure.

In view of the gap between patches, we create a buffer zone with a radius of 50 m for each patch, and it is assumed that two patches are adjacent when their buffer zones intersect. Hence, the adjacency relationship can be represented as the adjacency matrix below:

(2)$A=\left(\begin{array}{ccc}{a}_{11}& \cdots & a_{1n}\\ ⋮& \ddots & ⋮\\ a_{m1}& \cdots & a_{mn}\end{array}\right)$

where $a_{mn}$ is a parameter, and it is assumed that $a_{mn}=1$ if patches m and n are adjacent; otherwise, $a_{mn}=0$ if patches m and n are not adjacent.

3.1.3. Multiscale Features for Characterizing LUPs

This section defines and calculates features to characterize LUPs (Figure 5a). There are two challenges to be addressed: determining the types of features and calculating them.

For one thing, this study considers a series of features from 4 aspects to characterize LUPs, including the proportions of each kind of land use, landscape indices, the distance to the city center, and the kernel density of hospitals, schools, corporations, and commercial residents (Table 2).

Distances to the city center and the kernel density of various POIs can reflect the distribution of the population, GDP, and LST to a certain degree, while the proportions of diverse land uses and landscape indices can characterize land-use patterns from structural and functional perspectives, respectively [22]. Generally, some land-use surfaces, including artificial surfaces such as building roofs and paved roads, can aggravate the urban heat island effect, and the distribution of buildings and urban land-use structures can affect the absorption of solar radiation and the formation of air flows, then influence the heating effect [24]. Other kinds of land-use surfaces, such as some types of open urban parks and water bodies, can mitigate the urban heat island effect [25]. And the cooling effect is influenced by the size [26] and shape of the park.

Subsequently, 9 landscape indices covering 3 types (shape, aggregation, and diversity) are selected to characterize LUPs:

• (a). Shape indices. Shape indices are derived from the perimeter and area of a patch and quantify the complexity of landscape patch shapes [27], including the mean shape index (SHAPE_MN).

• (b). Aggregation indices. Aggregation indices measure the extent of spatial clustering or the dispersion of patches within a landscape including the patch density (PD), landscape shape index (LSI), contagion index (CONTAG), percentage of like adjacency (PLADJ), landscape division (DIVISION), and effective mesh size (MESH).

• (c). Diversity indices. Diversity indices originate from the concept of species diversity in ecology [28] and assess the landscape structure to reflect the diversity of land-use types [29], including Simpson’s Diversity Index (SIDI) and Simpson’s Evenness Index (SIEI). Thus, there is a lot of relevance between the proportions of land-use landscape indices and the LST. Among these features, only the proportion of park, industrial, commercial, and third-level residential land and landscape indices are variable, while the other features remain constant.

Furthermore, this section divides the land-use pattern into several of the same tiles and extracts these features tile by tile. Specifically, the distance to the city center is the Euclidean distance. Six landscape indices can be calculated using software called FRAGSTATS 4.2 [30] using the uniform tiles method. In this study, a 1 × 1 km tile was used to calculate the landscape patterns of the study area.

3.2. Multi-Objective Simulation for Diverse LUPs

3.2.1. Random Forest Regression for Quantifying LUPs’ Influences on Population, GDP, and LST

This section aims to analyze LUPs’ influences on the population, GDP, and LST (Figure 5b), and the key is to analyze the correlation between LUP features and the population, GDP, and LST. Aiming at the nonlinear interactions between the LUPs and the objectives, this study selected random forest regression, which is combined with several independent decision trees [31] and has the advantages of high accuracy and suitability for nonlinear data [32,33].

Basically, random forest regression is a machine learning-integrated algorithm that establishes and assembles multiple unrelated regression trees by randomly selecting samples and features and then integrates the results of all regression trees to improve the accuracy of the model [31]. The process of random forest regression is as follows: First of all, each regression tree in the random forest selects a part of the samples from the training data randomly with replacement [34]. The samples not selected are called out of bag (OOB) and can be used to assess regression errors. Then, for each split within a tree, a random subset of LUP features is selected as split attributes. Finally, the results of all regression trees are averaged as the final regression result. As for model performance evaluation, OOB errors based on OOB samples or R-Square can be selected to assess the errors of random forest regression. What is more, random forest regression can evaluate the importance of land-use pattern features to measure the influence strengths on the population, GDP, and LST. This section uses built-in random forest importance to evaluate feature importance, such as Gini importance. This study employs 100 trees for each random forest model and sets their maximum depth as 10 to ensure both convergence and efficiency of regressions.

To sum up, random forest can model the correlation between LUPs and optimization objectives. In this section, the study area is divided into several tiles of the same size (1 km × 1 km). Then, the population, GDP, and LST in each tile are calculated, and for each tile, several features of LUPs are defined and extracted. Finally, taking the value of all features in all tiles as the independent variable and the population, GDP, and LST in the corresponding tiles as dependent variables, 3 random forest regressions can be established to basically model the correlation between LUPs and the population, GDP, and LST.

3.2.2. Effectiveness Verification of Multi-Objective Simulation of LUPs

The mathematic relationship between the objectives and LUPs is the basis and difficulty of the LUP optimization method. In addition to choosing suitable features to characterize LUPs, the accuracy and robustness of the regression model also determine the reliability of the optimization results. Thus, it is significant to evaluate the accuracy and verify the effectiveness of multi-objective simulation. This study uses R-square (R2, Equation (3)) and Root Mean Squared Error (RMSE, Equation (4)) to evaluate the regression models (Figure 5b). The value of R-Square ranges from 0 to 1 and is closer to 1 when the model fits better. The RMSE represents the root of the average sum of squared errors and is close to 0 when the model fits better.

(3)$R^2=1-{\displaystyle \frac{{\displaystyle \sum _{i=1}^n{(y_i-\widehat{y})}^2}}{{\displaystyle \sum _{i=1}^n{(y_i-\overline{y})}^2}}}$

(4)$RMSE=\sqrt{{\displaystyle \frac{1}{n}}{\displaystyle \sum _{i=1}^n{(y_i-{\widehat{y}}_i)}^2}}$

where $y_i$ is the true value of the object, and ${\widehat{y}}_i$ is the predicted value of the object.

3.3. Spatiotemporal Constraints of LUP Optimization

3.3.1. Spatial Area Constraints of Land Uses

In order to coordinate the development of urban and rural areas and the intensive, economical, scientific, and rational use of land resources, planning standards of development land should be formulated to restrict the area of various types of land. An excessive area of a certain land-use category leads to the simplification of LUP functions. Thus, the area of each category should be constrained to a certain scope (Figure 5d), and the constraints are indicated by the formulas below.

(5)$\underset{¯}{p} S\le {\displaystyle \sum _{i\in I}{A}_{i,park}}\le \overline{p} S$

(6)$\underset{¯}{r} S\le {\displaystyle \sum _{i\in I}{A}_{i,\mathit{res}}}\le \overline{r} S$

(7)$\underset{¯}{i} S\le {\displaystyle \sum _{i\in I}{A}_{i,\mathit{indu}}}\le \overline{i} S$

where S is the total area of the study area; I is the collection of all land-use scenarios; $\underset{¯}{p}$ is the lower limit of the proportion of park area in the study area; $\underset{¯}{r}$ is the lower limit of the proportion of residential district area in the study area; $\underset{¯}{i}$ is the lower limit of the proportion of industrial land area in the study area; $\overline{p}$ is the upper limit of the proportion of park area in the study area; $\overline{r}$ is the upper limit of the proportion of residential district area in the study area; $\overline{i}$ is the upper limit of the proportion of industrial land area in the study area; $A_{i,park}$ is the area of parks in land-use scenario i; $A_{i,res}$ is the area of residential districts in land-use scenario i; and $A_{i,indu}$ is the area of industrial lands in land-use scenario i.

According to the code for the classification of urban land-use and planning standards of development land of China, which was released in 2011, construction land structure is shown in Table 3.

3.3.2. Neighborhood Dependency of Land Uses

As demonstrated in Section 3.4.2, this study selected a graph structure to present LUPs and represented the neighborhood relationships between patches by using an adjacent matrix. Then, how does this study use these neighborhood relationships during the process of optimization? How can neighborhood dependency be embodied?

In the traditional genetic algorithm, a crossover operator exchanges several random genes between two random chromosomes, allowing the genes to exploit particularly beneficial portions within a search space. And a mutation operator turns some random genes into new genes, exploring new regions within the search space [35]. Common crossover operators are single-point crossover, K-point crossover, uniform crossover, and polynomial crossover. However, according to Tobler’s First Law of Geography, “Everything is related to everything else, but near things are more related than distant things”. If a certain patch adjusts, the relationships with its neighboring patches change accordingly, resulting in the neighboring patches being more likely to renew. According to the representation methods of LUPs above, patches represented by neighboring genes on a chromosome may be far apart from each other. Aiming at this law, this study comes up with a modified crossover operator and mutation operator which take adjacency into consideration (Figure 5d). The new crossover operator exchanges the genes of the crossover site and its neighboring sites between two chromosomes. As shown in Figure 7a, assuming that patch 9 is selected as the crossover site, which is contiguous to patches 4–8, 10, and 11, we exchange the genes of patches 9–11 between two chromosomes.

While the traditional mutation operator only alters the gene of a randomly chosen mutation site, such as flip bit mutation, swap mutation, and inversion mutation, the new mutation operator also alters the gene of its neighboring sites. As shown in Figure 7b, we change the genes of patches 9–11.

3.3.3. Temporal Conversion Rules

The land-use optimization problem is a practical question, so this study should be combined with facts. In fact, some conversions between some land-use categories are hardly possible, such as the conversion from park to shanty town, which is unacceptable and may make the objects worse and then slow down the process of evolution. In order to guarantee the validity of the conversion, this study formulated type conversion rules (Figure 5d) by using a 6 × 6 matrix as follows:

(8)$rule=\left[\begin{array}{ccc}{r}_{11}& \cdots & r_{16}\\ ⋮& \ddots & ⋮\\ r_{61}& \cdots & r_{66}\end{array}\right]$

where $r_{mn}$ is a parameter, and it is assumed that $r_{mn}=1$ if the transformation from type m to n is allowable; otherwise, $r_{mn}=0$ if the transformation is not acceptable. The details of the temporal conversion rules are shown in Table 4.

3.4. Patch-Level and Neighborhood-Dependency Spatial Optimization Method (PNO)

3.4.1. Definition of Optimization Objective

On the one hand, as a city grows, changes in urban land cover and human activities in urban environments make urban areas warmer than the surrounding non-urbanized areas [36]. This phenomenon is the urban heat island and is the most obvious climate manifestation of urbanization [37]. On the other hand, Beijing city has absorbed more and more people who have settled here. Both of these effects have led to a heavy burden for the city. Thus, this study selects three optimization objectives including social–economic and ecological factors. The first objective is the population. In cities, some kinds of parcels can carry few people, so we can improve the population carrying capacity and make full use of the land through the optimization of the overall arrangement of land use. The second objective is the LST. In recent years, with the acceleration of the urbanization process, the relationship between humans and the environment has become increasingly tense, resulting in some eco-environmental problems and global warming. These problems have affected the living suitability seriously. Thus, aiming at the high land temperature, it is urgent to optimize the allocation and structure of land-use patterns to cool the land. Moreover, the LST is obtained from remotely sensed thermal infrared imagery, while the air temperature is only derived from a few meteorological stations. So, the LST is more continuous and can cover the entire cityscape [38,39]. The final objective is social economy. Economics development is still a non-negligible issue, and socio-economic objectives are always conflicting with ecological objectives, which will make some objectives worse while making other objectives better. Thus, this study made some improvements on the basis of the NSGA-II and came up with a patch-level and neighborhood-dependency spatial optimization method to solve the contradiction between economic and environmental factors.

3.4.2. Land-Use Optimization Methods Considering Patch Neighborhood Relationships

In order to solve the paradoxical problems above, this study came up with a patch-level and neighborhood-dependency spatial optimization method based on the genetic algorithm. The genetic algorithm is an intelligent optimization algorithm proposed by John H. Holland [40] which simulates the crossover operation and mutation operation in natural selection, evaluates the advantages and disadvantages of individuals by using the fitness function, and selects the individuals with higher fitness according to the principle of the survival of the fittest.

In order to sort all individuals by balancing all objectives, Srinivas and Deb [41] put forward a new algorithm based on the genetic algorithm named the Non-dominated Sorting Genetic Algorithm (NSGA). On the basis of this, Deb et al. [13] improved the NSGA and proposed NSGA-II. Elite strategy was introduced into this algorithm: parent population P_t and child population Q_t are merged into a population R_t, and non-dominated sorting is performed to divide the population into corresponding fronts. Then, the front with the highest non-dominated level is selected and put into the new-generation population. The pseudocode of the fast non-dominated sorting algorithm and the process of elite strategy are respectively shown in Algorithm 1 and Figure 8, while the NSGA-II procedure is shown in Figure 5c.

Algorithm 1 Fast non-dominated sorting (P)

for each p∈ P do

Sp = ∅           Used to store the members dominated by p     np = 0

for each q P do       if (p ≺ q) then        If p dominates q        Sp = Sp ∪ {q}        Add q to the set of solutions dominated by p       else if (q ≺ p) then        np = np + 1       Increase the domination counter of p     end     if np = 0 then        p belongs to the first front      prank = 1      F1 = F1 ∪ {p}  end  i = 1  while Fi ≠ ∅     Q = ∅           Used to store the members of the next front     for each p ∈ Fi do       for each q ∈ Sp do         nq = nq − 1         if nq = 0 then       q belongs to the next front        qrank = i + 1        Q = Q ∪ {q}       end     end     i = i + 1     Fi = Q  end

During the process of the genetic algorithm, the most important parts are the operations of crossover, mutation, and selection. And a patch’s land-use changing results in a change in local land-use patterns around it, making the patches more likely to alter. For example, the third-level residential districts around a newly opened mall would be more likely to be transformed to other land uses. Thus, this study uses the modified crossover and mutation operations which can embody the renewal and change law of adjacent patches.

Finally, taking the area within the fifth ring road of Beijing as the experimental area, setting the population size to 20, the number of iterations to 5000, and the probability of crossover and mutation to 0.8 and 0.1 separately, the experiment on land-use pattern optimization was conducted in five groups. Group A used the algorithm considering neighborhood relationships and temporal conversion rules. Group B used the algorithm only considering temporal conversion rules (without regard to neighborhood relationships). Group C used the algorithm only considering neighborhood relationships (without regard to temporal conversion rules). Group D used the algorithm without regard to neighborhood relationships and temporal conversion rules. Group E optimized the urban land-use pattern at the pixel level.

3.5. Evaluation of PNO Results

In Section 3.1, Section 3.2, Section 3.3 and Section 3.4, a patch-level and neighborhood-dependency spatial optimization method has been established, resulting in an optimized LUP that balances various objectives. In order to validate the effectiveness and scientific rigor of this optimized LUP, this section compares the optimized results with the LUP of the study area in 2022 (Figure 5e). An evaluation system has been established to measure effectiveness by examining the consistency between the optimized results and the 2022 LUP of the study area, focusing on both overall consistency and the consistency of changes.

Initially, overall consistency is measured by the proportion of identical land-use categories in corresponding patches between the optimized result and the 2022 LUP, which can be calculated by ArcGIS 10.8 overlap analysis. The consistency of changes is assessed by comparing the changes from the original 2020 LUP to both the optimized LUP and the 2022 LUP.

4. Results

After running the algorithm and performing data processing, we obtained results from three aspects, including the optimization results of LUP, the effectiveness evaluation results of PNO, and the advancement of the algorithm (Figure 9).

4.1. Optimization Results of LUPs

4.1.1. Quantitative Measurement for the Optimization Results

After running the algorithm with 5000 as the number of iterations, the optimal values of each objective were obtained and are shown in Table 5. In addition, taking the number of iterations as the x axis and the optimal GDP, population, and LST in the corresponding generation as the y axis separately, we can draw three curves of the objectives. In group A, the population, GDP, and LST witness growths of 4.72%, 2.75%, and −0.32% after optimization.

After optimization, we obtain 20 land-use patterns with different values for the population, GDP, and LST. Selecting the pattern with the lowest LST as the result of this experiment, the area statistics results of each land use before and after optimization and the area change rate are shown in Figure 10a, which illustrates that the area of parks, first-level residential districts, and second-level residential districts has increased after optimization, while areas of industrial lands, commercial lands, and third-level residential districts have been reduced. First-level residential districts have the highest growth rate (746.99%), followed by parks (2.92%), and the area of second-level residential districts has the lowest growth rate. Then, third-level residential districts have the largest decline (−53.43%), followed by industrial lands (−2.21%) and commercial lands (−2.25%). Qualitatively and quantitatively analyzing the change in land uses, the result is depicted and presented by the following Sankey Diagram (Figure 10b).

4.1.2. Analysis of Local Change in Land Uses

After analyzing the overall optimization, we focus on local changes in land uses. The original remote sensing image of the study area is presented in Figure 11a, the original and optimized LUPs are shown in Figure 11b,c, and the change from the original LUP to the optimized LUP is shown in Figure 11d. After comparison and analysis, noticeable changes are mainly concentrated in five zones.

To begin with, the first zone was mostly within the second loop which has the Palace Museum, and many national government departments and Hutongs, which is located in the most archaic and centered area of Beijing with a high population density and LST. Compared with the current land-use pattern, the third-level residential districts (L3-res) which are located to the north of the Palace Museum and a large area of institutions (Figure 11i) change into first-level residential districts (L1-res) and second-level residential districts (L2-res). Then, some second-level residential districts in the upper left corner of the zone i change into first-level residential districts and parks.

The second zone (Figure 11ii) is around the fourth loop, and the zone includes considerable areas of industries and Fengtai railway station and is divided by rails. Compared with the current land-use pattern, a considerable area of industries is converted into second-level residential districts.

The third zone (Figure 11iii) occurs in the northwest of the study area, named Zhongguancun. This zone has a large number of universities and commercial office buildings, so there are a lot of college students and commuters here. Compared with the current land-use pattern, some commercial lands are changed into first-level residential districts, and some second-level residential districts are changed into parks and first-level residential districts around the Zhongguancun Campus of the Chinese University of China.

The fourth zone (Figure 11iv) occurs to the south of the Palace Museum. Compared with the current land-use pattern, some third-level residential districts are converted to parks and first-level residential districts.

The final zone (Figure 11v) is to the west of Beijing Olympic Park, where many universities are, such as China Agricultural University and the China University of Geosciences. Compared with the current land-use pattern, some second-level residential districts near the Peking University Health Science Center Campus are changed to commercial lands, and some second-level residential districts near Beihang University and the University of Science and Technology Beijing are converted to first-level residential districts.

Combined with the locations and surroundings of the patches that need to be optimized and analyzed in the optimization, this study gives these reasonable explanations for these optimizations: (1) First-level residential districts are mainly low-rise buildings with a sparse population and lower impervious fraction; second-level residential districts are the typical high-density residential communities in Beijing, including high-rise buildings with a high-density population [42] and complete transportation, public service facilities, and more greening and bodies of water, absorbing heat from air during the day; and third-level residential districts are always correlated with low-quality and poor environments. The optimizations from third-level residential districts to second-level residential districts can effectively alleviate the situation of high LST, and the optimizations from second/third-level residential districts to first-level residential districts can help to increase populations and reduce the LST. (2) Urban parks consist of extensive green vegetation which can absorb the heat in the air by transpiration and block the sunlight from reaching the surface, as well as a considerable amount of water surfaces which can form a relatively cooler environment, along with the surrounding area within a certain range, and this cooling effect is quite significant and has a certain mitigation impact on the urban heat environment [43]. Parks also have a lower artificial surface rate than any other construction land which can reduce heat reflection. In addition, green spaces including parks provide residents with places for leisure activities and improve the quality of the residents’ lives. (3) Industrial lands can provide considerable jobs, create economic benefits, and then improve the economy of the surrounding patches. However, the proportion of industrial areas has a consistent and important positive correlation with the LST, demonstrating the contribution of industrial lands to the urban heat island effect [44]. So, appropriately reducing the area of industrial lands and increasing the area of vegetation in urban areas can reduce heat production from the source effectively.

4.1.3. Optimization Contributions to Sustainable Development

In the LST preferred pattern, the population rises from 14,129,368 to 14,796,691, increasing by 4.72% overall. Analyzing the population distribution at the 1 km tile scale specifically, the population increased in most of the tiles. For example, the population generally rose in the south and southwest of the study area (Figure 12i), which are thinly populated, and in some densely populated areas, as shown in Figure 12ii,iii, the population decreased.

As for the GDP, it rises from USD 388.56 billion to USD 399.25 billion, witnessing a growth of 2.75% as a whole after optimization. Analyzing the distribution of the GDP at the 1 km tile scale specifically, the GDP in the center of Beijing city and west of the city center (Figure 13i) was cut, while the GDP in the east, south, and north of the city center (Figure 13ii–iv) and in the south and north of the study area edge (Figure 13v,vi) saw a general increase. After optimization, the value of the GDP grew, and its distribution became scattered instead of being concentrated in the city center in the original pattern.

Finally, we turn to the LST. From an overall point of view, the average LST in the study area dropped by 0.12 degrees centigrade. Analyzing the distribution of the LST at the 1 km tile scale specifically, the cooling effect can achieve 8 to 9 degrees centigrade in the city center (Figure 14i), 2 to 8 degrees centigrade in the southwest of the city center (Figure 14ii) and in the east of the study area (Figure 14iii), and 1 to 2 degrees centigrade in the north of the city center (Figure 14iv). On the contrary, the LST increased by 0.72 to 7 degrees centigrade in the southeast (Figure 14v) and northwest (Figure 14vi) of the study area.

4.2. Results of the Effectiveness Evaluation

4.2.1. Multi-Objective Simulation Results Based on LUPs

In Section 3.2, the experiment trained a random forest regression model to analyze the LUPs’ influence on the population, GDP, and LST and then validated the effectiveness of multi-objective simulation by calculating the R-square (R²) and Root Mean Squared Error (RMSE). The evaluation results are listed in Table 6.

As shown in Table 6, the R-square values of the regression between land-use pattern features and the population, GDP, and LST are all more than 0.5, which means that these three regressions are all effective. Among them, the R-square value of regression between LUP features and LST is highest, and the RMSE value of regression between LUP features and LST is lowest, which means that the regression result between LUP features and LST is most accurate.

4.2.2. Effectiveness Evaluation of PNO Results

Comparing the proportion of the same land use of the corresponding patches between the optimization result and LUP of the study area in 2022, 43.8% of the area of the optimization result is consistent with the LUP in 2022, while it is 41.4% in the experiment without considering neighborhood relationships and temporal conversion rules. And compared with the changes from 2020 to 2022, there mainly exist three congruous areas consistent in the PNO optimization, occupying 25% of the optimized area while only accounting for 10.6% in the experiment without considering neighborhood relationships and temporal conversion rules. These results demonstrate that the optimized LUP generated by PNO is more in line with the current situation of land-use development.

As for the local zone, the first zone (Figure 15i) is located between the fourth and fifth ring roads, around Fengtai railway station, and contains a large number of factories and second-level residential districts. After optimization, some industrial lands in the zone (Figure 15iv) are converted to second-level residential districts, which is consistent with the land-use changes from 2020 to 2022.

The second zone (Figure 15ii) is located in the northeastern part of the study area called Wangjing, beside the fourth ring road, and contains abundant Internet enterprises and people. So, the land uses around there are mainly commercial lands, residential districts, and then parks with little areas of industrial land. Following optimization, the small amount of industrial land tagged in Figure 15v is transformed to parks and commercial land, which was also changed to parks and commercial land in the 2022 LUP. The increase in parks and decrease in industrial land both contribute to a reduction in the LST.

The last zone (Figure 15iii) is located in the western part of the study area, between the fourth and fifth ring roads, mainly composed of some residential districts and industrial lands, including cultural and creative parks, building materials companies, and concrete batching plants. After optimization, two industrial patches are turned into a second-level residential district and a park, which is partly consistent with the land-use changes from 2020 to 2022. The conversion from industrial land to a second-level residential district can not only increase the population in this patch but can also enhance the environment of the surrounding residential districts, then increase the occupancy of the surrounding residential districts.

4.3. Results of PNO Advancement Analysis

This section compares the results between group A and B–E to evaluate the effectiveness of PNO, which introduces two key improvements: the consideration of neighborhood relationships between land-use patches and the temporal conversion rules. The experiments are designed to analyze the individual and combined contributions of these improvements, demonstrating their impact on the overall performance of the optimization process and the advancement of PNO.

4.3.1. Without Considering Neighborhood Relationships

This section aims at conducting the experiment without considering neighborhood relationships and comparing the results with Section 4.1. In this section, the optimal population, GDP, and LST witness growths of 1.99%, 0.98%, and −0.12% after optimization, while this is 4.72%, 2.75%, and −0.32%, respectively, in the experimental group, which demonstrates that the performance of the genetic algorithm has improved each objective by 137.2%, 180.6%, and 166.7% in contrast to the experiment without considering neighborhood relationships formulated in this research.

Selecting the LST preferred pattern as the result of this comparative experiment, the original remote sensing image of the study area is presented in Figure 16a, the original and optimized LUPs of PNO are shown in Figure 16b,c, and the comparison result without considering neighborhood relationships is shown in Figure 16d. After comparison and analysis, noticeable differences are mainly concentrated in four zones.

To begin with, the first zone (Figure 16i) is mostly within the second loop which has the Palace Museum, and many national government departments and Hutongs, which is located in the most archaic and centered area of Beijing with a high population density and LST. Compared with the current land-use pattern, the third-level residential districts (L3-res) which are located to the north of the Palace Museum and a large area of institutions change into first-level residential districts (L1-res) and second-level residential districts (L2-res) in the experiment. Then, some second-level residential districts in the upper left corner of the zone i change into first-level residential districts and parks. But in the comparative experiment, the second- and third-level residential districts remain unchanged. A third-level residential district to the east of the Palace Museum was changed to a second-level residential district in the experiment but to commercial land in the comparative experiment.

The second zone (Figure 16ii) occurs in the northwest of the study area, named Zhongguancun, with a large number of universities and commercial office buildings, so there are a lot of college students and commuters here. Compared with the current land-use pattern, some commercial lands are changed into first-level residential districts in the experiment and are changed to parks in this comparative experiment.

The third zone (Figure 16iii) occurs to the south of the Palace Museum. Compared with the current land-use pattern, some third-level residential districts are changed to parks and first-level residential districts in the experiment and remain unchanged in this comparative experiment.

The last zone (Figure 16iv) is around the fourth loop, and the zone includes considerable areas of industries and Fengtai railway station and is divided by rails. Compared with the current land-use pattern, a considerable area of industries is converted into a second-level residential district in the experiment and remains unchanged in this comparative experiment.

In contrast to the comparative experiment without considering neighborhood relationships, the PNO method influences several contiguous patches within a certain range, instead of only some discrete patches. This is because the urban LUP is an integral and subtle system, where each kind of patch has its scope of influence. For example, industrial land generally has a negative impact on its surrounding districts, while parks have a positive influence. Thus, a patch’s land use altering may have a significant impact on its ambient local area. This demonstrates that the land use of adjacent patches generally changes in the result LUP of the PNO method.

4.3.2. Without Considering Temporal Conversion Rules

This section aims at conducting the experiment without considering temporal conversion rules and comparing the results with Section 4.1. In this section, the optimal population, GDP, and LST witness growths of 6.3%, 3.7%, and −0.64% after optimization, while this is 4.72%, 2.75%, and −0.32%, respectively, in the experimental group, which demonstrates that the performance of the genetic algorithm has decreased each objective by 33.47%, 34.55%, and 100% in contrast to the experiment without considering temporal conversion rules formulated in this research.

Selecting the LST preferred pattern as the result of this comparative experiment, the original remote sensing image of the study area is presented in Figure 17a, the original and optimized LUPs of PNO are shown in Figure 17b,c, and the comparison result without considering temporal conversion rules is shown in Figure 17d. After comparison and analysis, a noticeable difference is observed in the east of the study area, beside the east fifth loop. Compared with the current land-use pattern, a large area of second-level residential districts is changed to parks, and some parks are changed to third-level residential districts. Third-level residential districts are often high-density residential districts lacking public, transportation, and public service facilities, with low-quality and poor environments, which should be eliminated instead of being generated during the LUP optimization. Accordingly, this conversion is invalid and completely divorced from reality.

4.3.3. Without Considering Both Neighborhood Relationships and Temporal Conversion Rules

This section aims at conducting the experiment without considering temporal conversion rules and neighborhood relationships and comparing the results with Section 4.1. In this section, the optimal population, GDP, and LST witness growths of 0.59%, 0.77%, and −0.08% after optimization, while this is 4.72%, 2.75%, and −0.32%, respectively, in the experimental group, which demonstrates that the performance of the genetic algorithm has improved each objective by 700%, 257.1%, and 300% in contrast to the experiment without considering temporal conversion rules and neighborhood relationships formulated in this research.

Selecting the LST preferred pattern as the result of this comparative experiment, the original remote sensing image of the study area is presented in Figure 18a, the original and optimized LUPs of PNO are shown in Figure 18b,c, and the comparison result without considering both neighborhood relationships and temporal conversion rules is shown in Figure 18d. After comparison and analysis, noticeable differences are mainly concentrated in three zones.

To begin with, the first zone (Figure 18i) was mostly within the second loop which has the Palace Museum, and many national government departments and Hutongs, which is located in the most archaic and centered area of Beijing with a high population density and LST. Compared with the current land-use pattern, a large area of third-level residential districts around the Palace Museum is changed to second-level residential districts, and a little part is converted to first-level residential districts in the experiment. However, in this comparative experiment, the third-level residential districts remain the same.

The second zone (Figure 18ii) is in the southwest corner of the study area where considerable areas of industrial lands are changed to second-level residential districts, and a second-level residential district near the railway is changed to a park in this comparative experiment, but in the experiment this stayed unchanged.

The third zone (Figure 18iii) occurs to the south of the Palace Museum. Compared with the current land-use pattern, some third-level residential districts are changed to parks in the experiment. However, in this comparative experiment without considering neighborhood relationships and temporal conversion rules, these third-level residential districts stay unchanged. What is more, several patches are changed to third-level residential districts, which is invalid and completely divorced from reality.

4.3.4. Pixel-Level vs. Patch-Level Optimization

Like many previous studies on pixel-level LUP optimization, there are various categories of fragments not next to each other in the result of this comparative experiment Figure 19, dividing the complete previous patches of land use into different parts, reducing spatial connectivity, spatial relationships, and connectivity between patches and aggravating the issue of land-use fragmentation. Meanwhile, the optimized results at pixel level are more difficult to integrate with actual management measures, weakening the feasibility and effectiveness of policy implementation.

5. Discussion

5.1. How Does PNO Resolve the Three Issues in LUP Optimization?

As demonstrated in Section 1, there are three hard issues in LUP optimization: (1) the urban land-use system is an integrated system, and the objectives are mutually conditioned, so the pattern may change a lot even if a patch changed; (2) the mathematic relationship between objectives and LUPs is unclear, and the most suitable method to analyze the relationship between objectives and LUPs is unknown, which makes it hard to calculate the objectives; and (3) current optimization methods optimize on a pixel-by-pixel basis and ignore the spatial dependencies between the neighboring patches. And PNO resolves these three issues with the following methods.

Initially, as it is widely used in current research about LUP optimization, NSGA-II is also used in the PNO to search for the optimal LUPs, which was brought forward by Deb et al. [13] and based on GA and NSGA. Using an efficient non-dominated sorting approach with a computational complexity of O(MN²), NSGA-II significantly reduces the computational burden compared to earlier algorithms. Then, NSGA-II uses an elitist approach by combining the parent and offspring populations and selecting the best solutions, ensuring that the best solutions are preserved across generations, improving convergence towards the Pareto-optimal front. In addition, NSGA-II employs a crowding distance mechanism to maintain diversity in the population. To sum up, NSGA-II can handle the hard problems and balance the contradictions of each objective and has been extensively studied and tested in LUP optimization studies.

Subsequently, this study selects twenty-eight features from four aspects, including the proportion of land use, distance to city center, kernel density of POIs, and landscape index, then divides the study area into the same tiles (1 km × 1 km), calculates these features, and makes zonal statistics on the population, GDP, and LST in each tile. Finally, taking features as independent variables and objectives as dependent variables, this study embraces random forest to analyze LUP features’ influence on the population, GDP, and LST. Through effectiveness verification, it is proved that the multi-objective simulation using random forest is effective and feasible.

Ultimately, this study represents the land-use patches as a one-dimensional array ordered by ID and a graph structure, where the former acts as a chromosome for evolution in the GA, and the latter represents the neighborhood relationships of patches. Then, new crossover and mutation operators which work both on the random patch and its adjacent patches are put forward in Section 3.4.2.

5.2. Pros and Cons of PNO in Urban Land-Use Planning

Previous studies of land-use pattern optimization focused on the pixel scale and ignored the neighborhood relationships of land uses, which however had three issues: (1) the pixel-by-pixel optimization process may destroy the integrity of the patch and cause land-use fragmentation which may decrease biodiversity and affect the integrality of urban land uses and activities, restricting the application of land optimization models in land and space planning, assistance, and decision making; (2) the neighborhood relationships of land uses, which can reflect internal correlation, were not taken into consideration during the optimization process; and (3) differently, this study proposes a patch-level and neighborhood-dependency spatial optimization method; thus, this study contributes an integrated LUP optimization model based on the genetic algorithm which can optimize LUPs at patch level and consider neighborhood dependency among patches. The pros and cons of this study are detailed as follows.

First and foremost, patches are the basic units constituting an LUP instead of pixels during land planning and land optimization. Government and planning departments usually manage and implement policies at the patch level. Patch-level optimization results are more easily integrated with actual management measures, improving the feasibility and effectiveness of policy implementation [45,46,47]. The pixel-by-pixel optimization process may destroy the integrity of the patch and cause land-use fragmentation. For example, a scenario where a residential area emerges inside a commercial patch (e.g., a shopping mall) after optimization is out of touch with reality. Meanwhile, since patch-level models deal with larger units of data, the data requirements and computational complexity are significantly lower compared to pixel-level models. This study focuses on patches, optimizing the LUP while preventing land-use fragmentation and providing more feasible and effective support for developing policies and decision making.

Additionally, this study takes the spatial interaction and neighborhood effects of patches into consideration during the optimization process. An LUP is a complex integrated system and always exhibits spatial auto-correlation and neighborhood effects for the spatial interactions between land-use categories itself [47]. When the use of one patch changes (e.g., from agricultural to commercial use), it can trigger changes in economic activities in the surrounding patches. A new commercial center might attract more commercial and residential development, leading to changes in the land-use types of the surrounding patches. This method takes this effect into consideration and implements it by modifying the crossover and mutation operations.

In summary, this work plays a significant role in advancing research on urban land use and could aid in scientific land-use optimization and planning. Although this study shows strong performance, there still exist several limitations. For example, some historical buildings carrying profound historical significance and humanistic value, such as Hutong, are real cultural heritage, and should be protected and repaired rather than torn down, but they are merely classified as third-level residential districts together with other old residential buildings and converted into another land use. In response to this issue, more reasonable and scientific constraints should be further considered for these historical buildings combined with historical protection building data in future studies. For example, these historic buildings could be set as restricted optimization areas, where urban land-use optimization is not allowed. Then, this study optimizes LUPs exclusively through transforming the land-use categories of patches rather than their location, shape, or size, but this always happens in urban planning. The above results in low consistency in the change in the optimization results relative to the original LUPs and actual change. Thus, a complex algorithm should be proposed.

6. Conclusions

Previous research on LUP optimization focused on the pixel level and seldom considered the impact on neighboring patches and the surrounding local land-use pattern if a certain patch’s land use changes. However, government and planning departments usually manage and implement policies at the patch level. The mismatch between the research on LUP optimization and practical applications led to a decrease in the feasibility of the optimization method; differently, this study represents LUPs at patch level and LUPs’ adjacency by using a graph data structure. Furthermore, this study proposes a patch-level and neighborhood-dependency spatial optimization method, which is the first work to consider neighborhood relationships in the process of optimization and modify the crossover and mutation operators in the genetic algorithm instead of optimizing LUPs patch by patch. Accordingly, this study contributes to improving the feasibility and effectiveness of LUP optimization.

Two conclusions can be drawn from this study. First of all, this work presents a patch-level and neighborhood-dependency spatial optimization method, which can make trade-offs between socio-economic development and sustainability objectives to optimize the LUP of the study area. Overall, this method performs better than previous spatial optimization methods. From the perspective of target distribution, areas with high population density have experienced a decrease in population density and population transfer to areas with a low population. And areas with a high LST have experienced a decrease in the LST.

Subsequently, this study evaluates the effectiveness of PNO results by comparing PNO results with the land-use classification in 2022. In reference to the 2022 land-use map, it is found that 43.8% of the study area in the PNO results matches the actual land-use classification in 2022, which is an improvement over the previous optimization results, where only 41.4% of the area matches. Then, this study also conducts a detailed comparison between the land-use changes in the PNO results and the actual land-use changes that occurred from 2020 to 2022. Compared with the actual land-use changes, 25% of the optimized area in the PNO results is aligned with the actual land-use changes from 2020 to 2022, exceeding previous spatial optimization results using NSGA-II with a consistency of 10.6%. Thus, the PNO is more effective and applicable in land-use planning and urban renewal.

In summary, these conclusions are completely data-driven, and the patch-level and neighborhood-dependency spatial optimization method can assist in urban land-use planning, decision making, and making trade-offs between socio-economic development and sustainability. While there exist serval limitations in this work, they will be further investigated in our future work.

Author Contributions

Conceptualization, Yuhan Cheng and Xiuyuan Zhang; data curation, Xiaoyan Dong; formal analysis, Yuhan Cheng; funding acquisition, Xiuyuan Zhang and Shihong Du; investigation, Xiaoyan Dong; methodology, Yuhan Cheng and Xiuyuan Zhang; project administration, Xiuyuan Zhang, Qi Zhou, and Shihong Du; resources, Xiuyuan Zhang and Xiaoyan Dong; software, Yuhan Cheng; supervision, Xiuyuan Zhang, Qi Zhou, and Shihong Du; validation, Yuhan Cheng and Xiuyuan Zhang; visualization, Yuhan Cheng and Xiaoyan Dong; writing—original draft, Yuhan Cheng; writing—review and editing, Xiuyuan Zhang. All authors have read and agreed to the published version of the manuscript.

Data Availability Statement

The data presented in this study are openly available on FigShare at https://doi.org/10.6084/m9.figshare.27175473 accessed on 12 January 2025.

Conflicts of Interest

The authors declare no conflicts of interest.

Footnotes

Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Figures and Tables

Figure 1. Differences from existing optimization methods, including (a) current pixel-by-pixel optimization method; (b) current patch-by-patch optimization method; and (c) optimization method proposed in this study.

Figure 2. Distribution of districts and the study area in Beijing (within 5th ring road).

Figure 3. Experimental data, including (a) land-use data in 2020, (b) land-use data in 2022, (c) road dataset overlapped with the satellite image, (d) population, (e) GDP, and (f) LST dataset within 5th loop of Beijing.

Figure 4. Point-of-interest (POI) data within 5th loop of Beijing which are sorted into 4 categories.

Figure 5. Technology roadmap, including 5 parts: (a) representation method of LUPs: describes the method of representing LUPs and their neighborhood relationships using an array and a graph structure, followed by the definition and calculation of features from different perspectives to characterize LUPs; (b) multi-objective simulation for diverse LUPs: explains how to quantify the influences of LUPs on population, GDP, and LST using random forest regression, along with verifying the effectiveness of multi-objective simulation for LUPs; (c) patch-level and neighborhood-dependency spatial optimization method (PNO): the main body of the improved optimization algorithm; (d) spatiotemporal constraints of LUP optimization: defines several constraints of LUP optimization, including spatial area, neighborhood dependency, and temporal conversion rules; and (e) effectiveness evaluation of PNO results: details the evaluation process for the effectiveness of PNO.

Figure 6. Graph structure for presenting LUPs, including (a) land use of local area composed of patches A–J and (b) its graph structure.

Figure 7. Crossover operator (a) and mutation operator (b) considering patch neighborhood relationships. Number 1–13 represent several patches. The red one represents the crossover or mutation site and the blue ones represent its adjacent patches.

Figure 8. Process of elite strategy. In the process, Pt and Qt are respectively parent population and child population in generation t. They are merged into a population Rt. and F1–F6 represent the fronts with the non-dominated level 1–6. Pt+1 is new population in next generation.

Figure 9. The results from three aspects, including the optimization results of LUP, the results of the effectiveness evaluation, and the results of PNO advancement analysis.

Figure 10. (a) Area changes in land uses and (b) chord diagram which describes the changes in land uses.

Figure 11. Result of experimental group A, including (a) remote sensing image overlapped with five frames which represent five selected example areas for visual analysis of local land use change in the five dashed frames below; (b) original land uses; (c) optimized result; (d) change areas; (e–h) visual analysis of local land use change in area i; (i–l) visual analysis of local land use change in area ii; (m–p) visual analysis of local land use change in area iii; (q–t) visual analysis of local land use change in area iv; and (u–x) visual analysis of local land use change in area v. Each set of visual analysis of local land use change contains local remote sensing image, original land uses, optimized result, and change areas.

Figure 12. Population distribution before (a) and after (b) the optimization and changes (c). Frames i–iii are selected example areas for detailed discussion.

Figure 13. GDP distribution before (a) and after (b) the optimization and changes (c). Frames i–vi are selected example areas for detailed discussion.

Figure 14. LST distribution before (a) and after (b) the optimization and changes (c). Frames i–vi are selected example areas for detailed discussion.

Figure 15. Effectiveness evaluation results of PNO results by comparing optimization result and LUP in 2022, including three local areas i–iii. (a–d) effectiveness evaluation results of PNO results in example area i; (e–h) effectiveness evaluation results of PNO results in example area ii; and (i–l) effec-tiveness evaluation results of PNO results in example area iii. Each set of effectiveness evaluation results contains remote sensing image, original land uses, optimization land uses, and land uses in 2022. Frames iv–vi are the local areas where the optimization is the same as the actual changes.

Figure 16. Comparison between the optimization result of experimental groups A and B, including (a) remote sensing image overlapped with four frames which represent four selected example areas for visu-ally comparing the optimizing results in the four dashed frames below.; (b) original land uses; (c) optimized result of PNO; (d) result without considering neighborhood relationships; (e–h) com-parison of optimizing results in example area i; (i–l) comparison of optimizing results in example area ii; (m–p) comparison of optimizing results in example area iii; and (q–t) comparison of op-timizing results in example area iv. Each set of comparisons contains local remote sensing image, original land uses, optimized result of PNO, and result without considering neighborhood rela-tionships.

Figure 17. Comparison between the optimization result of experimental groups A and C, in-cluding (a) remote sensing image overlapped with one frame which represents an example area for visually comparing the optimizing results in the dashed frame below.); (b) original land uses; (c) optimized result of PNO; (d) result without considering temporal conversion rules; (e–h) com-parison of optimizing results (local remote sensing image, original land uses, optimized result, and result without considering temporal conversion rules) in example area i.

Figure 18. Comparison between the optimization result of experimental groups A and D, including (a) remote sensing image overlapped with three frames which represent three selected example areas for visually comparing the optimizing results in the three dashed frames below.; (b) original land us-es; (c) optimized result of PNO; (d) result without considering both neighborhood relationships and temporal conversion rules; (e–h) comparison of optimizing results in example area i; (i–l) comparison of optimizing results in example area ii; and (m–p) comparison of optimizing results in example area iii. Each set of comparisons contains local remote sensing image, original land uses, optimized result of PNO, and result without considering neighborhood relationships and temporal conversion rules.

Figure 19. Comparative results by using patch-level and pixel-level method, including (a) remote sensing image; (b) original land uses; (c) optimized result at patch level; and (d) optimized result at pixel level.

Table 1

Land-use category for each code.

Code	Land-Use Category	Descriptions
1	Forest	Forests and grassland
2	Water	Rivers, sea, lakes, ponds, etc.
3	Unused lands	Undeveloped wastelands
4	Transportation	Rural roads, streets, highways, railways, railway stations, airports, and parking lots
5	Parks	Land that is used for health, sports, and other facilities
6	Industrial lands	Industrial, mining, and storage lands
7	Institutions	Science, education, culture, and health areas
8	Commercial lands	Commercial offices, restaurants, hotels, shopping malls, shops, supermarkets, etc.
9	L1-res	High-end and low-density residential land with complete public, transportation, and public service facilities and good environment
10	L2-res	High-density residential land with relatively complete public, transportation, and public service facilities and good environment
11	L3-res	Dilapidated housing, shanty towns, temporary housing, etc., lacking public, transportation, and public service facilities and with poor environment
12	Farmland	Land that is used for farming

Table 2

Features characterizing local land-use patterns.

Type	Name	Descriptions
Proportion	[Proportion]	Proportions of each land-use category in each slice.
Distance	[Euclidean_Dis]	Euclidean distance from palace museum to each slice.
Kernel density	[Kernel_Den]	Kernel density of four types of POIs, including hospitals, schools, corporations, and commercial residents.
Landscape indices	[PD]	Patch density measures the number of land-use patches per unit area.
	[LSI]	Measures the landscape shape index.
	[SHAPE_MN]	Shape indices are based on the perimeter and area of a patch and represent the complexity of landscape patch shapes.
	[CONTAG]	The contagion index indicates the degree of aggregation of various landscape types and their spatial distribution characteristics within the landscape.
	[PLADJ]	Proportion of like adjacencies measures the percentage of adjacent landscape types in the total adjacent landscape.
	[DIVISION]	Division measures the degree of landscape fragmentation.
	[MESH]	Effective mesh size measures the areas of diverse land uses per unit area.
	[SIDI]	Simpson’s Diversity Index measures the richness of land-use categories in the landscape.
	[SIEI]	Simpson’s Evenness Index measures the evenness of land-use categories’ frequency distribution.

Table 3

Planned construction land structure.

Category Name	Proportion of Urban Construction Land
Residential land	25.0–40.0
Administration and public services land	5.0–8.0
Industrial land	15.0–30.0
Transportation facilities land	10.0–30.0
Green space	10.0–15.0

Table 4

Temporal conversion rules of land-use categories.

Land-Use Category	Park	Industrial Land	Commercial Land	L1-Res	L2-Res	L3-Res
Park	0	1	1	1	1	0
Industrial land	1	0	1	1	1	0
Commercial land	1	0	0	1	1	0
L1-res	1	1	1	0	1	0
L2-res	1	1	1	1	0	0
L3-res	1	1	1	1	1	0

Table 5

Result of objectives in experimental group.

Objective	Value		Change (%)
	Before Optimization	After Optimization
GDP	14,129,368.34	14,796,691.25	4.72
Population	388.56	399.25	2.75
LST	36.86185731	36.74368767	−0.32

Table 6

R-Square and MSE of each regression model.

Objective	R2 (%)	RMSE
GDP	0.60	0.96
Population	0.51	802.2
LST	0.66	1.008

References

1. Wang, W.; Jiao, L.; Jia, Q.; Liu, J.; Mao, W.; Xu, Z.; Li, W. Land use optimization modelling with ecological priority perspective for large-scale spatial planning. Sustain. Cities Soc.; 2021; 65, 102575. [DOI: https://dx.doi.org/10.1016/j.scs.2020.102575]

2. Zhang, T.; Sun, L.; Yao, L.; Rong, J. Impact Analysis of Land Use on Traffic Congestion Using Real-Time Traffic and POI. J. Adv. Transp.; 2017; 2017, 7164790. [DOI: https://dx.doi.org/10.1155/2017/7164790]

3. Qin, P.; Wang, L. Job Opportunities, Institutions, and the Jobs-Housing Spatial Relationship: Case Study of Beijing. Transp. Policy; 2017; 81, pp. 331-339. [DOI: https://dx.doi.org/10.1016/j.tranpol.2017.08.003]

4. Cao, K.; Zhang, W.; Wang, T. Spatio-temporal land use multi-objective optimization: A case study in Central China. Trans. GIS; 2019; 23, pp. 726-744. [DOI: https://dx.doi.org/10.1111/tgis.12535]

5. Du, S.; Xiong, Z.; Wang, Y.-C.; Guo, L. Quantifying the multilevel effects of landscape composition and configuration on land surface temperature. Remote Sens. Environ.; 2016; 178, pp. 84-92. [DOI: https://dx.doi.org/10.1016/j.rse.2016.02.063]

6. Sahebgharani, A. MULTI-OBJECTIVE LAND USE OPTIMIZATION THROUGH PARALLEL PARTICLE SWARM ALGORITHM: CASE STUDY BABOLDASHT DISTRICT OF ISFAHAN, IRAN. J. Urban Environ. Eng.; 2016; 10, pp. 42-49. [DOI: https://dx.doi.org/10.4090/juee.2016.v10n1.042049]

7. Sicuaio, T.; Zhao, P.; Pilesjo, P.; Shindyapin, A.; Mansourian, A. Sustainable and Resilient Land Use Planning: A Multi-Objective Optimization Approach. ISPRS Int. J. Geo-Inf.; 2024; 13, 99. [DOI: https://dx.doi.org/10.3390/ijgi13030099]

8. Schlager, K.J. A LAND USE PLAN DESIGN MODEL. J. Am. Inst. Plan.; 1965; 31, pp. 103-111. [DOI: https://dx.doi.org/10.1080/01944366508978151]

9. Rahman, M.M.; Szabó, G. Multi-objective urban land use optimization using spatial data: A systematic review. Sustain. Cities Soc.; 2021; 74, 103214. [DOI: https://dx.doi.org/10.1016/j.scs.2021.103214]

10. Wang, G.; Han, Q.; de Vries, B. The multi-objective spatial optimization of urban land use based on low-carbon city planning. Ecol. Indic.; 2021; 125, 107540. [DOI: https://dx.doi.org/10.1016/j.ecolind.2021.107540]

11. Cao, K.; Liu, M.; Wang, S.; Liu, M.; Zhang, W.; Meng, Q.; Huang, B. Spatial Multi-Objective Land Use Optimization toward Livability Based on Boundary-Based Genetic Algorithm: A Case Study in Singapore. ISPRS Int. J. Geo-Inf.; 2020; 9, 40. [DOI: https://dx.doi.org/10.3390/ijgi9010040]

12. Li, X.; Parrott, L. An improved Genetic Algorithm for spatial optimization of multi-objective and multi-site land use allocation. Comput. Environ. Urban Syst.; 2016; 59, pp. 184-194. [DOI: https://dx.doi.org/10.1016/j.compenvurbsys.2016.07.002]

13. Deb, K.; Pratap, A.; Agarwal, S.; Meyarivan, T. A fast and elitist multiobjective genetic algorithm: NSGA-II. IEEE Trans. Evol. Comput.; 2002; 6, pp. 182-197. [DOI: https://dx.doi.org/10.1109/4235.996017]

14. Balling Richard, J.; Taber John, T.; Brown Michael, R.; Day, K. Multiobjective Urban Planning Using Genetic Algorithm. J. Urban Plan. Dev.; 1999; 125, pp. 86-99. [DOI: https://dx.doi.org/10.1061/(ASCE)0733-9488(1999)125:2(86)]

15. Lu, Y.; Wang, S.; Zhao, Y.; Yan, C. Renewable energy system optimization of low/zero energy buildings using single-objective and multi-objective optimization methods. Energy Build.; 2015; 89, pp. 61-75. [DOI: https://dx.doi.org/10.1016/j.enbuild.2014.12.032]

16. Wang, A.; Xiao, Y.; Liu, C.; Zhao, Y.; Sun, S. Multi-objective optimization of building energy consumption and thermal comfort based on SVR-NSGA-II. Case Stud. Therm. Eng.; 2024; 63, 105368. [DOI: https://dx.doi.org/10.1016/j.csite.2024.105368]

17. Thakur, K.; Maity, S.; Nielsen, P.; Pal, T.; Maiti, M. A 3D multiobjective multi-item eco-routing problem for refrigerated fresh products delivery using NSGA-II with hybrid chromosome. Comput. Ind. Eng.; 2024; 198, 110644. [DOI: https://dx.doi.org/10.1016/j.cie.2024.110644]

18. Liu, H.; Yan, F.; Tian, H. Towards low-carbon cities: Patch-based multi-objective optimization of land use allocation using an improved non-dominated sorting genetic algorithm-II. Ecol. Indic.; 2022; 134, 108455. [DOI: https://dx.doi.org/10.1016/j.ecolind.2021.108455]

19. Luan, C.; Liu, R.; Zhang, Q.; Sun, J.; Liu, J. Multi-objective land use optimization based on integrated NSGA–II–PLUS model: Comprehensive consideration of economic development and ecosystem services value enhancement. J. Clean. Prod.; 2024; 434, 140306. [DOI: https://dx.doi.org/10.1016/j.jclepro.2023.140306]

20. Wilcox, B.A.; Murphy, D.D. Conservation Strategy: The Effects of Fragmentation on Extinction. Am. Nat.; 1985; 125, pp. 879-887. [DOI: https://dx.doi.org/10.1086/284386]

21. Rotem-Mindali, O. Retail fragmentation vs. urban livability: Applying ecological methods in urban geography research. Appl. Geogr.; 2012; 35, pp. 292-299. [DOI: https://dx.doi.org/10.1016/j.apgeog.2012.08.003]

22. Zhang, X.; Du, S.; Du, S.; Liu, B. How do land-use patterns influence residential environment quality? A multiscale geographic survey in Beijing. Remote Sens. Environ.; 2020; 249, 112014. [DOI: https://dx.doi.org/10.1016/j.rse.2020.112014]

23. Zhao, N.; Liu, Y.; Cao, G.; Samson, E.L.; Zhang, J. Forecasting China’s GDP at the pixel level using nighttime lights time series and population images. GISci. Remote Sens.; 2017; 54, pp. 407-425. [DOI: https://dx.doi.org/10.1080/15481603.2016.1276705]

24. Gago, E.J.; Roldan, J.; Pacheco-Torres, R.; Ordóñez, J. The city and urban heat islands: A review of strategies to mitigate adverse effects. Renew. Sustain. Energy Rev.; 2013; 25, pp. 749-758. [DOI: https://dx.doi.org/10.1016/j.rser.2013.05.057]

25. Taha, H. Urban climates and heat islands: Albedo, evapotranspiration, and anthropogenic heat. Energy Build.; 1997; 25, pp. 99-103. [DOI: https://dx.doi.org/10.1016/S0378-7788(96)00999-1]

26. Cao, X.; Onishi, A.; Chen, J.; Imura, H. Quantifying the cool island intensity of urban parks using ASTER and IKONOS data. Landsc. Urban Plan.; 2010; 96, pp. 224-231. [DOI: https://dx.doi.org/10.1016/j.landurbplan.2010.03.008]

27. Wu, J.; Liu, Z.M. Effect of habitat fragmentation on biodiversity: A review. Chin. J. Ecol.; 2014; 33, pp. 1946-1952.

28. Roberts, F.S. Measurement of Biodiversity: Richness and Evenness. Mathematics of Planet Earth: Protecting Our Planet, Learning from the Past, Safeguarding for the Future; Kaper, H.G.; Roberts, F.S. Springer International Publishing: Cham, Switzerland, 2019; pp. 203-224.

29. Feng, Y.; Liu, Y.; Tong, X. Spatiotemporal variation of landscape patterns and their spatial determinants in Shanghai, China. Ecol. Indic.; 2018; 87, pp. 22-32. [DOI: https://dx.doi.org/10.1016/j.ecolind.2017.12.034]

30. McGarigal, K.S.; Cushman, S.; Neel, M.; Ene, E. FRAGSTATS: Spatial Pattern Analysis Program for Categorical Maps; University of Massachusetts: Amherst, MA, USA, 2002.

31. Breiman, L. Random Forests. Mach. Learn.; 2001; 45, pp. 5-32. [DOI: https://dx.doi.org/10.1023/A:1010933404324]

32. Ouma, Y.O.; Nkwae, B.; Odirile, P.; Moalafhi, D.B.; Anderson, G.; Parida, B.; Qi, J. Land-Use Change Prediction in Dam Catchment Using Logistic Regression-CA, ANN-CA and Random Forest Regression and Implications for Sustainable Land–Water Nexus. Sustainability; 2024; 16, 1699. [DOI: https://dx.doi.org/10.3390/su16041699]

33. Roy, B. A machine learning approach to monitoring and forecasting spatio-temporal dynamics of land cover in Cox’s Bazar district, Bangladesh from 2001 to 2019. Environ. Chall.; 2021; 5, 100237. [DOI: https://dx.doi.org/10.1016/j.envc.2021.100237]

34. Breiman, L. Bagging predictors. Mach. Learn.; 1996; 24, pp. 123-140. [DOI: https://dx.doi.org/10.1007/BF00058655]

35. Datta, D.; Fonseca, C.M.; Deb, K. A multi-objective evolutionary algorithm to exploit the similarities of resource allocation problems. J. Sched.; 2008; 11, pp. 405-419. [DOI: https://dx.doi.org/10.1007/s10951-008-0073-9]

36. Voogt, J.A.; Oke, T.R. Thermal remote sensing of urban climates. Remote Sens. Environ.; 2003; 86, pp. 370-384. [DOI: https://dx.doi.org/10.1016/S0034-4257(03)00079-8]

37. Fan, H.; Sailor, D.J. Modeling the impacts of anthropogenic heating on the urban climate of Philadelphia: A comparison of implementations in two PBL schemes. Atmos. Environ.; 2005; 39, pp. 73-84. [DOI: https://dx.doi.org/10.1016/j.atmosenv.2004.09.031]

38. Klok, L.; Zwart, S.; Verhagen, H.; Mauri, E. The surface heat island of Rotterdam and its relationship with urban surface characteristics. Resour. Conserv. Recycl.; 2012; 64, pp. 23-29. [DOI: https://dx.doi.org/10.1016/j.resconrec.2012.01.009]

39. Nichol, J.E.; Fung, W.Y.; Lam, K.-s.; Wong, M.S. Urban heat island diagnosis using ASTER satellite images and ‘in situ’ air temperature. Atmos. Res.; 2009; 94, pp. 276-284. [DOI: https://dx.doi.org/10.1016/j.atmosres.2009.06.011]

40. Holland, J.H. Adaptation in Natural and Artificial Systems: An Introductory Analysis with Applications to Biology, Control and Artificial Intelligence; MIT Press: Cambridge, MA, USA, 1992.

41. Srinivas, N.; Deb, K. Muiltiobjective optimization using nondominated sorting in genetic algorithms. Evol. Comput.; 1994; 2, pp. 221-248. [DOI: https://dx.doi.org/10.1162/evco.1994.2.3.221]

42. Yao, L.; Xu, Y.; Zhang, B. Effect of urban function and landscape structure on the urban heat island phenomenon in Beijing, China. Landsc. Ecol. Eng.; 2019; 15, pp. 379-390. [DOI: https://dx.doi.org/10.1007/s11355-019-00388-5]

43. Xu, L.H. A study on thermal environment effect of urban park landscape. Acta Ecol. Sin.; 2008; 28, pp. 1702-1710.

44. Lu, Y.; Yue, W.; Huang, Y. Effects of Land Use on Land Surface Temperature: A Case Study of Wuhan, China. Int. J. Environ. Res. Public Health; 2021; 18, 9987. [DOI: https://dx.doi.org/10.3390/ijerph18199987] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/34639286]

45. Liang, X.; Guan, Q.; Clarke, K.C.; Liu, S.; Wang, B.; Yao, Y. Understanding the drivers of sustainable land expansion using a patch-generating land use simulation (PLUS) model: A case study in Wuhan, China. Comput. Environ. Urban Syst.; 2021; 85, 101569. [DOI: https://dx.doi.org/10.1016/j.compenvurbsys.2020.101569]

46. Parker, D.C.; Manson, S.M.; Janssen, M.A.; Hoffmann, M.J.; Deadman, P. Multi-Agent Systems for the Simulation of Land-Use and Land-Cover Change: A Review. Ann. Assoc. Am. Geogr.; 2003; 93, pp. 314-337. [DOI: https://dx.doi.org/10.1111/1467-8306.9302004]

47. Verburg, P.H.; Schot, P.P.; Dijst, M.J.; Veldkamp, A. Land use change modelling: Current practice and research priorities. GeoJournal; 2004; 61, pp. 309-324. [DOI: https://dx.doi.org/10.1007/s10708-004-4946-y]