Practical Efficient Regional Land-Use Planning Using Constrained Multi-Objective Genetic Algorithm Optimization

Abstract

Practical efficient regional land-use planning requires planners to balance competing uses, regional policies, spatial compatibilities, and priorities across the social, economic, and ecological domains. Genetic algorithm optimization has progressed complex planning, but challenges remain in developing practical alternatives to random initialization, genetic mutations, and to pragmatically balance competing objectives. To meet these practical needs, we developed a Land use Intensity-restricted Multi-objective Spatial Optimization (LIr-MSO) model with more realistic patch size initialization, novel mutation, elite strategies, and objectives balanced via nominalizations and weightings. We tested the model for Dapeng, China where experiments compared comprehensive fitness (across conversion cost, Gross Domestic Product (GDP), ecosystem services value, compactness, and conflict degree) with three contrast experiments, in which changes were separately made in the initialization and mutation. The comprehensive model gave superior fitness compared to the contrast experiments. Iterations progressed rapidly to near-optimality, but final convergence involved much slower parent–offspring mutations. Tradeoffs between conversion cost and compactness were strongest, and conflict degree improved in part as an emergent property of the spatial social connectedness built into our algorithm. Observations of rapid iteration to near-optimality with our model can facilitate interactive simulations, not possible with current models, involving land-use planners and regional managers.

1. Introduction

Global population growth and rapid unbalanced urban expansion in the 21st century [1] has caused many social–environmental issues such as forest degradation, soil erosion, biodiversity decline, and irreversible changes in some land uses [2]. In addition, pressures from population growth have increased land-use demand for transport, industries, commercial buildings, and housing, resulting in conflicts between different land uses within restricted geographic zones [3]. Effective urban planning involves the suitable allocation and resolution of conflicting land uses that must also balance regional land resource demand and supply [4]. Therefore, systematic land-use optimization involves balancing allocation tradeoffs across complex land-use objectives that must also align with regional development constraints and aspirations; preferably, it also involves planners interacting with models to explore the impacts of uncertainties and alternate formulations. Therefore, our aim with this study is to detail improvements in land-use planning methodology with a realistic application demonstration that is much better resolved than previous studies.

Current genetic algorithm (GA)-based land-use optimization models are generally focused on balancing competing objectives and evolving changes across potentially viable alternate patterns to enhance the reliability of the optimization [5,6,7]. However, practical implementation challenges remain in producing final land-use optimizations that meet the requirements of constraints and priorities of regional policies. First, initializations of cells with random land use results in highly fragmented patches, unsuitable land-use patterns, and computational inefficiencies in the labored iterations toward a potentially pre-mature and non-optimal solution [8,9,10,11]. Second, some internal operators in GA (e.g., crossover and mutation) crucially affect the compatibility of alternative land uses and the model’s efficiency [5,12,13,14]; thus, further consideration is needed on the rationality of land-use conversion required to resolve patterns. Thirdly, the balance between conflicting land uses requires careful scaling and weighting to avoid an unbalanced comparison of competing objectives with differing value profiles [14,15,16,17,18]—for example, economic differences between very high value residential land and relatively low value forests.

Land-use optimization is essentially about reaching compromises in land-use co-existence, a modus vivendi, that in an agreed sense best balances regional community aspirations for society, ecology, and the economy [8]. However, differences among localities lead to various combinations of conflicting objectives. Therefore, tradeoffs are required, as it is impossible to completely meet the needs of all objectives everywhere. At present, planners’ objectives are mainly divided into three competing categories: economic benefit, environment benefit, and social benefit [7]. Economic benefit is mainly measured by the unit output or by the conversion cost of each land use [19,20], either by Gross Domestic Product (GDP) allocation or conversion cost, respectively. Within land-use models, economic factors can be evaluated for maximum housing capacity, maximum social equity, and minimum taxes [18,21,22,23]. Among various ecological factors, Ecosystem Services Value (ESV), namely the value of ecosystem service and natural capital estimated by economic rules, is widely used due to its clear scientific significance and ease of calculation [24,25,26]. In addition, Net Primary Productivity (NPP), soil erosion, and pollution control can be chosen if necessary to accord with environmental protection policies [27,28]. For the social benefit measurement, compatibility or incompatibility of land use is most often adopted [29,30,31], while compactness and conflict degree (between neighboring land use localities) are also used for social objectives [3,32,33].

Present land-use optimization algorithms are generally of three types: linear optimization, cellular automata (CA), and intelligent algorithm [33]. Linear optimization estimates the aspatial optimal area of each land use and cannot deal with spatial arrangements [9,34]. CA uses spatial-neighbor conversion rules to time step the changing land-use pattern, but it cannot be used to explore relationships among economic, ecological, and other factors [35]. The intelligent algorithm overcomes the shortcomings of above two algorithms by solving Non-deterministic Polynomial (NP) hard problems, which by definition are decision problems that are not solved deterministically step-by-step in time, but whose solutions, derived by non-deterministic methods (non-deterministic Turing machine), are easy to check from its fitness value. Although intelligent algorithms may not converge to optimal solutions in some cases, due to the complexity of the problem, near-optimal solutions are also meaningful [10]; hence, this is the approach adopted for this study.

Intelligent algorithms utilize some form of expert knowledge in decision making [36,37] and include the simulated annealing (SA) algorithm [19], ant colony optimization (ACO) method [38], particle swarm optimization (PSO) [39], and GA [40]. Although SA can avoid local optimizations traps, Aerts et al. (2005) demonstrated that GA performs better than SA for solution time and compactness [41]. The hybrid PSO requires setting a constant proportion of land use as a constraint before space allocation, whereas this is not required in GA [42]. The essence of ACO is to modify the transition probability of cells through positive feedback loops [43]. However, due to the conflicting nature of objectives, it is not clear which objective should be used as feedback during the evolutionary process. These general considerations led us to choose GA for our study in view of its global search capability and robustness.

Global optimization with GA algorithms iteratively improves (non-deterministically) the evolution of solutions by naturally selecting higher quality (“fitness”) genes, that is to say child solutions that are iteratively generated from parent solutions by sequential operations of selection, crossover, and mutation [3,25,44]. Iteration here means non-deterministic algorithmic iteration (not necessarily time-stepping) where it is possible to check if the “fitness” of the iterated solution is improving or not.

The reliability of the solution in avoiding premature near-optimal, or local optimum traps, can be enhanced by improved internal operators or by combining with other algorithms. For instance, to enhance the accuracy of the model, Cao et al. (2012) designed a three-step mutation operator based on patch-based diversity, boundary-based, and elimination of solutions to satisfy constraints [10]. Ma et al. (2018) used adaptive probability and elite strategy to avoid local optimization traps [32]. Huang et al. (2014) combined GA and CA to obtain the Pareto optimal solution for dynamically expanding cities [18]. Although those methods can improve model reliability, we propose that optimizations more suited for practical implementation can be obtained if minimum patch size constraint, a common requirement of land-use planning, is implemented in the initialization. For example, some initialization improvements have been made such as generating parent solutions with cell conversion probability [25]. The essential argument with initialization improvements is that the initial parent land-use patterns are closer to the optimal solution (than random allocation) and therefore should improve the search efficiency and final land-use distribution. In contrast, random initial seeds will form highly fragmented land uses, causing difficulties in achieving viable spatial agglomerations in new growth patches. Moreover, in crossover or mutation, cell changes mostly depend on the number of adjacent cells of the same type, while other factors (e.g., social and spatial benefits) are rarely considered to improve efficiency.

In solving the above problems, and in pursuit of an intensive land-use pattern with less isolated cells and enhanced aggregation ability, we developed a model named ‘Land use Intensity-restricted Multi-objective Spatial Optimization’ (LIr-MSO), which improves initialization, mutation, and elite strategies across multi-objectives. We selected two economic objectives, one ecological objective and two social objectives, while government policies and planning standards are implicit in specific constraints. Our study provides three distinct contributions: the first is the adjustment of the initialization strategy to improve final optimization reliability by constraining minimum patch sizes for different land uses; the second is an improvement of the mutation and elite strategy to make it more efficient; lastly, in the application, we demonstrate how to scale nominalizations to balance out competing objectives, and the effect of implicit tradeoffs in the model.

2. Materials and Methods

2.1. Study Area and Data Resource

Dapeng Community, located in the east of Shenzhen city, was utilized as the demonstration study area. In length, it is about 14.2 km north to south and 18.8 km east to west, with a total coverage area of about 88.3 km². More than 45,000 people live and work in this area (up to 2018). To facilitate the optimization computation, the whole region is rasterized into 1878 × 1418 = 2,663,004 cells, with each cell being 10 × 10 m in size (Figure 1). This number of cells is a substantial increase over past models and an indication of our progress to more realistic models than past simulations; for example, 16,779 cells were used by Cao, Huang, Wang, and Lin [10] and 1600 cells were used in the early model by Stewart, Janssen, and Herwijnen [32].

Each cell is given a unique land-use type from the following Classes (uppercase character abbreviation): Residential land (R), Commercial land (C), Industrial land (I), Arable land (Ar), Aquaculture land (Aq), Forest land (F), Other land (O), Beach (B), Public land (P), and Water (W).

Social and economic data were derived from Dapeng’s new 13th five-year plan and the statistical yearbook 2017 of Dapeng District [45]. Land-use types (construction land, forest land, arable land, aquaculture land, beach, and water) were visually interpreted and then digitized from remotely sensed images (GF-1, 2017) [46], and finally, construction land was further divided into commercial land, industrial land, public land, and residential land, by calibrating with points of interest from Baidu map [47]. Then, the land-use data were rasterized to the formulated spatial extent and resolution. Digital Elevation Model (DEM) data from SRTM (Shuttle Radar Topography Mission) [48] were clipped and interpolated to the rasterized land-use grid.

2.2. Model for Land-Use Optimization: LIr-MSO

As discussed, fitness impact factors/values are necessary to tradeoff patterns and to determine one (the fittest) that optimally meets the constraints and requirements of the land-use objectives. These factors are economic efficiency, ecological protection, spatial allocation, and regional development policies. According to the inherent properties of such factors, economic, ecological, and spatial factors are often regarded as objectives to be optimized, and social factors are implicitly represented by constraints and allowable juxtapositions of alternate land use in applying government policies and relevant social standards.

Figure 2 shows the workflow procedure of LIr-MSO, which involves initialization, fitness calculation, genetic evolution, and then feedback to the fitness calculation from the elite strategy until the land-use pattern is optimized or the maximum iteration is reached. First, initial patterns are obtained by setting minimum patch sizes to different land-use types. Second, the fitness of the pattern is calculated for the customized and nominalized objectives and constraints. Third, roulette is adopted to select candidates based on fitness. Then, if allowed by their crossover probabilities, any two patterns from candidates are selected for crossover to generate a new pattern. Next, the pattern is further optimized by using a two-step, block-based and point-based, spatial mutation operator, and a judgment is made to determine whether to go to the next step or to the next loop by comparing the fitness between the pattern and its parents. Finally, the elite strategy determines whether the fitness of the optimal pattern in this generation is better than that in the previous generation, and the final solution is obtained after customized iterations.

2.2.1. Objectives

Since land-use change involves both static and dynamic adjustments to intrinsic and conversion value respectively, a single index cannot comprehensively represent the economic benefits. Therefore, we choose both conversion cost (dynamic) and GDP (static) to represent economic benefits. For the environment benefit, ESV is chosen due to its wide use. Moreover, both compactness and conflict degree are adopted for evaluating social benefit because of their advantages in measuring the social compatibility of land use.

• Objective 1: Minimization of Conversion Cost

The conversion cost is the cost to redevelop the land-use type of one cell to another land use [49]. According to the literature and expert knowledge, a matrix of transfer coefficients is used to represents per-unit area conversion cost between all possible land-use conversions [20,28,50]. The sum of the product of transfer coefficient and changed area represents the total conversion cost for an existing land-use pattern. The transfer coefficient ranges from 0 to 1. The formula for objective 1 ($O\mathit{1})$ is as follows:

$$\mathrm{Min}O\mathit{1} = {{\displaystyle \mathit{\sum }}}_{s=1}^{ K}{{\displaystyle \mathit{\sum }}}_{t=1}^{ K}{B}_{\mathit{st}}{\mathit{TA}}_{\mathit{st}}.$$

(1)

B_st represents the transfer coefficient from existing land-use type s to type t; TA_st denotes the transfer area from land-use type s to type t; $K$ is the number of land-use types.

• Objective 2: Maximization of GDP

GDP varies according to the regional economic development level, and its real value involves regional, social, and environmental components partitioned amongst the different land uses. In this paper, the unit GDP value for each land use is calculated by using land-use area and GDP statistics. The total GDP for a given land-use pattern was obtained by summing the product of the unit GDP value and corresponding land-use area. The formula for objective 2 ($O\mathit{2})$ is as follows:

$$\mathrm{Max}O\mathit{2} = {{\displaystyle \mathit{\sum }}}_{k=1}^{ K}{U}_k{A}_k.$$

(2)

U_k is unit GDP value for land-use type k; A_k is total changed area of land-use type k; $K$ is the number of land-use types.

• Objective 3: Maximization of ESV

ESV is generally interpreted as the value of services provided by ecosystems and ecological processes utilized by humans [51]. The ecological benefit for objective 3 ($O\mathit{3})$ is estimated as follows (Costanza et al., 1997) [26]:

$$\mathrm{Max} O\mathit{3} = {{\displaystyle \sum }}_{k=1}^K{A}_k{V}_k.$$

(3)

V_k represents the ESV of the land-use type k per unit area; A_k represents the total changed area of land-use type k; $K$ is the number of land-use types.

• Objective 4: Maximization of Compactness

A reasonable pattern normally arranges each land-use cell in coherent companionship with similar and complementary uses to achieve social use objectives. Land-use compactness is considered a key factor that enhances living convenience and land utilization efficiency. Methods for measuring compactness include adjacency-based clustering, perimeter-area based compactness, and block aggregation [10,19,29]. After reviewing these approaches, we chose the popular perimeter-area based compactness for its simplicity in application. In this method, compactness is defined by the perimeter divided by the square root of area; higher values indicate less compact aggregation of the same land-use type and vice versa for lower values. The formula for objective 4 ($O\mathit{4})$ is as follows:

$$\mathrm{Min}O\mathit{4} = {{\displaystyle \mathit{\sum }}}_{k=1}^K\frac{{{\displaystyle \mathit{\sum }}}_{h=1}^{H_k}\frac{P_{hk}}{ \sqrt{A_{hk}}}}{H_k}.$$

(4)

P_hk is perimeter of the h-th patch in land-use type k; A_hk is area of the h-th patch in land-use type k; H_k is number of patches of land-use type k; $K$ is the number of land-use types.

• Objective 5: Minimization of Conflict Degree

To improve the suitability of the living environment, the incompatibility between different land-use types should be minimized. The incompatibility is often caused by conflicts from urban noise, pollution, land-use density, and other factors. Thus, a quantitative measure of conflict is necessary. We use degree of unsuitability between two adjacent land-use types as a measure of the degree of conflict (measured on a scale from 0 to 8), which should be low to ensure harmony between land uses. The formula for objective $O\mathit{5}$ is calculated using the eight neighbors (j), around each center cell (i), as illustrated below (Figure 3).

The formula for objective 5 (O5) is then calculated as follows:

$$\mathrm{Min} O\mathit{5} = {{\displaystyle \sum }}_{i=1}^{ N}{{\displaystyle \mathit{\sum }}}_{j=\mathit{1}}^{ 8}{{\displaystyle \sum }}_{k=1}^K{{\displaystyle \sum }}_{l=1}^K{X}_{ik}{X}_{jl}{C}_{kl}.$$

(5)

X_ik is a binary variable. X_ik = 1 when cell i is developed into land-use type k, otherwise X_ik = 0; C_kl is the conflict degree between two adjacent cells with land-use type k and l in the changed area; $K$ is the number of land-use types.

2.2.2. Constraints

Specific constraints that align with current planning policies and land-use standards are as follows:

1. A unique land-use type is assigned to each cell to retain cohesiveness and to avoid conflicts; this is implicit in the modified mutation algorithm;

2. The upper and lower area limits of each land use are constrained to comply with planning criteria and government policies, thereby minimizing high fragmentation;

3. Land-use change is excluded within the ecological red line delineated by the local government, to enhance sustainable living;

4. Terrain with slopes steeper than 25 degrees can only be used as forest land to prevent environmental impacts such as erosion.

2.2.3. Procedures of the LIr-MSO Model

To improve the model’s efficiency, some enhancements are made in the procedure, which comprises six steps: initialization, fitness calculation, selection, crossover, mutation, and elite strategy. The multi-objective spatial optimization problem can be defined as:

$$\mathrm{Max} \sum Oi, i\in [1, 5]$$

(6)

where minimum objectives (i = 1 to 5) are transformed into maximum objectives after normalization.

Step 1: Initialization

During initialization, patches within the changeable area are randomly generated to obtain N solutions to form an initial population. Each land-use alternative is regarded as being analogous to a chromosome in which genes correspond to grid cells, whose value represents a special land-use type. Thus, genetic interaction among different chromosomes can lead to different land-use patterns and optimization results. The commonly used initialization distributes changed cells in a random pattern consisting of many isolated cells. This fragmentation affects the efficiency of the algorithm, as many agglomeration trials are necessary across the whole region, leading to many iterations slowing progress toward the final solution [21].

As also as suggested by Maitland [52] in the minimal urban structure theory, the frequency of nodes represented by obvious urban forms (e.g., crossroads) constitutes the fundamental atomic structure of a city. Maitland [52] suggested that 200 m (4 ha in area for a 200 × 200 m block) could be used as the minimal unit structure of a city to guide larger development and construction. Thus, to obtain a representative and practical land-use pattern, the expansion of land use should be initiated with a minimal land-use-type unit as an atomic building block (with different atoms/block sizes for different land uses) for larger genetic land-use structures.

Guided by the minimal block scale theory of Maitland [52] and urban planning standards in China, the minimum area for residential land is set as 2 ha (comprising 200 cells, each of 100 m²), and for commercial land, industrial land, and public land, the minimum area is 0.5 ha (50 cells), respectively. As demonstrated in the following experiments, such initiation will help avoid fragmentation and the phenomenon of large patches engulfing small patches in the following operators (e.g., crossover and mutation) whereby some initially generated patches might be merged in the final solution, and also some isolated cells of impractical size might still exist.

We generated initial patterns depending on the conversion probability of each land use, the formula is as follows:

$$C{R}_k={{\displaystyle \sum }}_{t=1 }^{K}C{P}_{tk}{A}_{tk}.$$

(7)

CP_tk is the conversion probability from land-use type t to land-use type k; A_tk is the total changed area from land-use type t to type k; $K$ is the number of land-use types.

Step 2: Fitness Calculation

Currently, two methods are widely adopted: one is the Pareto-based method and the other is the weighted method. The former resolves optimal solution sets without compromising different objectives; hence, the relative priority of objectives does not determine optimality. The latter assigns relative weights, determined heuristically, that prioritize each objective so that the weighted sum of the objectives represents a single fitness measure that can be used as an index of optimality. We adopted the latter to calculate the fitness, as it allows planners and managers to tailor their priorities through the selection of appropriate weights. To account for magnitude differences of objectives, we normalized each into a range from 0 to 1. Considering the inherent characteristics of objectives, we adopted three normalization strategies: Objectives 2 (GDP) and 3 (ESV) need to be maximized and scaled according to the normalization formula:

$$Z\left(V\right) = \frac{{ V - V}_{\mathit{min} }}{V_{\mathit{max}} {- V}_{\mathit{min}}}.$$

(8)

Objectives 1 (conversion cost) and 5 (conflict degree) need to be minimized so their nominalization is as follows:

$$Z\left(V\right) = \frac{V_{\max }- V}{V_{\mathit{max}} {- V}_{\mathit{min}}}.$$

(9)

V is the raw objective fitness; $V_{\mathit{max}}, V_{\mathit{min}}$ are the respective maximum and minimum fitness for that objective.

For objective 4 (compactness), the normalization formula allows for compactness having a fixed minimum value (for a square block = 4 * R/$\sqrt{R^2}$ = 4, where R is the length of the block):

$$N\left(X\right)= \frac{X_{min}}{X}.$$

(10)

X is the raw objective fitness; $X_{\mathit{min}}$ is the minimum fitness nominally 4 for a square.

Then, the weighted fitness is obtained as follows:

$$\mathit{fit}\left(x\right){ = w}_1{\mathrm{Z}\left(\mathrm{O}1\right) + w}_2{Z\left(\mathrm{O}2\right) + w}_3{Z\left(O3\right) + w}_5{\mathrm{Z}\left(\mathrm{O}4\right) + w}_{4}N\left(\mathrm{O}5\right)$$

(11)

where fit(x) represent the fitness of the xth chromosome; w_i represents the weight of objective i. The sum of weights is 1.

Step 3: Selection

We use roulette in the selection operation so that the probability of each chromosome being selected is directly proportional to its fitness. The higher the fitness of the chromosome, the greater the probability of its selection in the next generation. The probability formula for cell i (1 to N cells) is as follows:

$${P(x}_i{ ) = \mathit{fit}(x}_i) / ({{\displaystyle \sum }}_{i=1 }^N{\mathit{fit}(x}_i\left)\right).$$

(12)

The following algorithm details the roulette selection of x_i:

R = rand ()

(13)

$$s = \left\{\begin{array}{cc}1\hfill & {R<P(x}_i) \hfill \\ 0\hfill & {R>P(x}_i)\hfill \end{array}\right..$$

(14)

R is a random number between 0 and 1. If the selection score (s) is equal to 1, x_i will be admitted to participate in the next step; otherwise, it will enter into the next selection to determine s; where the other possibility is that it will not be selected (s = 0).

Step 4: Crossover

The purpose of the crossover is to reduce the fragmentation of land-use patterns [25]. The procedure for crossover [7,33] is as follows: after calculating the crossover probabilities of parents, parent 1 (main parent) is selected along with parent 2 (subsidiary parent); a portion of cells from the parents are randomly selected according to their crossover probabilities (determined as frequencies of occurrence). For example, in Figure 4, supposing two cells, one from parent 1 and another from parent 2, have different land-use types (a and b respectively) in the same position; a 3*3 window around the selected cell is used to separately count the number of cells in parent 1 that belong to a and b. The selected cell in parent 1 (a) will be replaced by b if the new number of b cells will be more than a after replacement.

Step 5: Mutation

Mutation promotes solution diversity by making changes to genes in the chromosome. To avoid local optimum traps and promote the compactness of land use, we use a two-step spatial mutation operator:

(1) Block-based Mutation (BM) promote compactness as follows: after customizing BM’s probability and selecting a parent, a portion of cells in the parent will be randomly selected according to the BM probability. Then, a 3*3 window around the selected cell is used to count the number of cells of each land-use type. Next, the conversion cost from one land-use type to another is calculated according to the transfer coefficient. Finally, cells within the window are replaced with the type that has the lowest conversion cost.

(2) Point-based Mutation (PM) is performed after BM to reduce geographic dispersion (Figure 5). The main procedure is as follows: after customizing PM’s probability and selecting a parent, a portion of cells in the parent will be randomly selected according to PM’s probability. Suppose the land-use type of the selected cell is a, a 3*3 window around the selected cell is used to change the cell if the number of a cells is less than other types; for example, if type b is the most, the selected cell will be replaced by b. An alternative way to think of this mutation is that of a median filter in which the center cell is replaced with the most frequently occurring within the window. Viewed in this way, the mutation operation acts as a spatial smoothing filter.

Step 6: Elite Strategy

Elite strategy, or “survival of the fittest”, is often adopted to promote convergence by retaining the fittest individuals in each generation [52]. However, iterations in the elite strategy can converge slowly toward sub-optimal fitness amongst a juvenile (pre-mature) set of parents. In order to avoid this local trap and to promote global convergence toward maturity and optimality, a replacement mode is adopted (Figure 6). In this mode, the best pattern is conserved as an additional resource to be used when the current generation fitness is less than that in the previous generation. Assume that A(n) is the parent population, A(n + 1) is a newly generated population, and n is the number of the iteration. In the replacement mode, we discuss three cases: the best fitness in A(n) is higher than A(n + 1), the best fitness in A(n) is equal to A(n + 1), and the best fitness in A(n) is lower than A(n + 1). The crossover is adopted only if the maximum fitness in A(n + 1) is not greater than that in A(n), thus ensuring that the fitness keeps increasing.

Step 7: Terminate the algorithm

We ensure that the number of generated chromosomes is consistent with the initial number in the whole process. The algorithm will not achieve the optimal result until N-times’ iterations (500 in the case study) are finished or the fitness is deemed to reach convergence (fitness does not vary over many iterations).

3. Results

3.1. Objective Quantification and Constraints

Basic parameters for the objectives were determined from several papers, published reports, and expert knowledge, and they are given in the following tables.

Table 1. Transfer coefficient, from a scale of 0 (change not possible) to 1 (no restrictions in swapping over) for changing land-use types. Abbreviations are as follows for the land-use type: Residential land (R), Commercial land (C), Industrial land (I), Arable land (Ar), Aquaculture land (Aq), Forest land (F), Other land (O), Beach (B), Public land (P), and Water (W). Columns are the source land use and rows are the changed land use.

Change to\Change from	R	C	I	O	F	Ar	P
R	0	0.44	0.32	0.18	0.22	1	0.24
C	0.23	0	0.31	0.18	0.22	1	0.27
I	0.29	0.37	0	0.18	0.22	1	0.45
O	1	1	1	0	1	1	1
F	0.89	0.92	0.89	0.7	0	0.3	0.85
Ar	0.83	0.85	0.82	0.7	0.4	0	0.81
P	0.35	0.51	0.33	0.34	0.18	1	0

: Transfer coefficients between two land-use types were compiled from published studies [20,28,50] and presented in Table 1 where rows are the transferred land use, while the columns are the source land use. Water bodies, beaches, and aquaculture land are unchangeable areas; hence, their transfer coefficients are 0.

In

Table 2. GDP (Gross Domestic Product) and ESV (Ecological Service Value) per unit area for different land-use types. Land-use type abbreviations are as per Table 1.

Land-Use Type	GDP per Unit Area (RMB/m2)	ESV per Unit Area
R	30,467	0
C	699	0
I	4508	0
P	0	0
Ar	0	7.9
Aq	4.9	7.9
F	0.001	28.12
O	0	1.39
B	0	28.12
W	0	45.35

, the unit GDP of different land uses were determined from remote sensing interpretation and regional economic statistics. Table 2 also lists the ESV for each land use as determined from the current land use in China [51]. The GDP of residential land is used as the maximum for nominalization, whilst for ESV, the value for forest is used as water is not changeable, despite its ESV being higher.

In

Table 3. Conflict degree between different land uses. The matrix is symmetrical so that conflict degree does not depend on the horizontal directionality of the spatial juxtaposition. Land-use type abbreviations are as per Table 1.

	R	C	I	O	F	Ar	P
R	0	4	8	5	7	8	0
C	4	0	6	5	7	8	2
I	8	6	0	5	7	8	7
O	5	5	5	0	1	1	5
F	7	7	7	1	0	2	6
Ar	8	8	8	1	2	0	8
P	0	2	7	5	6	8	0

, the conflict degree for any two land uses is summarized from the present literature and expert knowledge [28,30]. Moreover, area constraints for each land use are shown in

Table 4. Regional area constraints for residential and industrial land use in hectares (ha).

Land Use Type	Lower Area Limits (ha)	Upper Area Limits (ha)
Residential land	No less than 410	No more than 660
Industrial land	No less than 250	No more than 500

.

The unchangeable area is shown in Figure 7a, including reserved villages (e.g., Dapeng Ancient City) and ecological preservation regions. Figure 7b shows a slope-constraint map that identifies potential construction lands where the slope is less than 25 degrees.

3.2. Implementation and Evaluation

Every objective uniquely affects corresponding benefits, making it difficult to determine which objective is the most important. In practice, expert knowledge of the tradeoffs from weightings is necessary, which is possibly in conjunction with scenario experimentation. As this is not the main aim of this study, we assume for demonstration purposes that regional development is balanced so that every objective has equal weighting; thus, their weights are all set 0.2 [8,24,36]. Furthermore, the probabilities of crossover, BM-based mutation, and PM-based mutation are 0.06, 0.05, and 0.3, respectively as determined from past studies [18,25,31]. A total of 100 chromosomes are initially set, and a maximum of 500 iterations (derived from preliminary experimentation) are performed.

One nominalization difficulty is with the GDP valuations which differ by over seven orders of magnitude (≈10⁷) between the highest value and lowest values (30,467 for Residential to 0.001 for Forest respectively). In preliminary experiments, we tested the nominalization of Equation (8) and found that the nominal GDP fitness value of the optimal pattern was not much more than 0.08 compared to, for example, compactness at over 0.9, leading to GDP playing a very minor role in the overall fitness contribution. With such a large disparity in GDP values, we resorted to taking log₁₀ of the numerator and the denominator of Equation (8) to raise the nominal value of GDP (to more than 0.8).

Before achieving the final solution, the summed value of objectives that contribute to fitness will increase until the fitness reaches a convergence (or the maximum designated iterations is reached). The change of fitness with iterations (Figure 8) shows an initial phase of rapid increase in fitness up to about 21 iterations and then a much slower phase of gradual fitness increase, including almost static phases where the fitness remains relatively stable through successive iterations. Shorter periods of these static phases are also noticeable in the rapid phase, so it is not a unique facet of the slow-change phase. The reason for the static phase is that, in this model, if the value in the parent generation is greater than that in the offspring, the optimal pattern in the parent generation is preserved. The curve appears stable after 226 iterations—more than 10 times the number of iterations required for the first rapid phase.

The maximum fitness is about 0.8520 after 226 iterations, while the fitness at the end of the rapid phase is 0.8422 after rising from 0.7645 in 21 iterations. So, the fitness change in the rapid phase, relative to the range from start to end, is about eight times greater than the slow phase. What is more, the rate of change (relative change per iteration) is about 77 times greater during the rapid phase, which is an indication that the LIr-MSO model is very efficient at achieving near-optimal fitness within what is a relatively short period of iterations compared to much longer iterations reported in other studies (e.g., Cao, Huang, Wang and Lin [10]) for much smaller models.

Five representative regions are selected to describe the differences between the initial and the optimal patterns (Figure 9). Compared to the initial pattern (

Table 5. Comparison of the change in area (km²) between initial pattern and optimal multi-objective optimal pattern for the different land uses. Abbreviations are described in Table 1.

Objs\Land Use (km2)	C	I	O	F	Ar
Initial	2.3894	6.5216	2.0319	67.2372	1.7036
LIr-MSO	4.9478	8.1120	1.8262	62.8576	1.6818
Changed (%)	107.07%	24.39%	−10.12%	−6.51%	−1.28%
Objs\Land Use (km2)	Aq	W	B	R	P
Initial	1.2413	1.3628	0.2270	4.8931	0.6966
LIr-MSO	1.2413	1.3628	0.2270	5.2754	0.7726
Changed (%)	-	-	-	7.81%	10.91%

, Figure 9 and Figure 10), the best solution has more commercial, industrial, public, and residential land, while forest, arable, and other land decrease. The main reasons for this change are the development needs of the economy. This suggests that there is still potential for further development in the study area while protecting the ecological environment. In the optimal pattern, residential land increased by 7.81%, mainly in Regions 2 and 3. The increase in residential areas near existing patches increased compactness. Commercial land increased by 107.07%, some patches merged in Region 2, and some newly generated patches appear in Region 1 and Region 5. Some increased commercial land along the coast suggests an increased effective utilization of coastal resources for regional GDP. Most of the reduced forest land was converted for further development to commercial, industrial, and residential land. In addition, in both the initial and optimal patterns, the industrial land situated in the west and east coasts is always surrounded by forest land, which can effectively reduce the impact of industrial pollution on human settlements. Moreover, industrial land use increased by 24.39%, mainly in Region 2 and Region 5, and several complementary commercial and residential patches appear within the industrial region, which can provide living space to workers to improve their quality of life.

In order to demonstrate the reliability of the model, three contrast experiments were conducted: (1) without initialization revision, namely patches generated without area constraint (condition 1); (2) without mutation revision (condition 2), namely all eight neighbors around the central cell converted to central cell’s land-use type; (3) without initialization and mutation revision (condition 3).

In Condition 1 (Figure 11a), although the full model has increased fitness by 1.74%, if no area constraint is adopted, land use will mainly change through patch expansion, contraction, and merging; thus, new patches will not be generated. In addition, this results in expanded, extended patches with a more rounded appearance that is most likely due to the compactness objective, as also evident in the comprehensive optimal pattern. While optimal from a compactness perspective, in practice, rounding off straight transport-related boundaries may be unacceptable in certain situations.

The pattern from Condition 2 (Figure 11b) suggests that the proposed mutation operator better handles isolated patches because in this simulation, there are some fragments, and the comprehensive fitness simulation is better by 2.05%.

In Condition 3 (Figure 11c), the fitness of the full model is better by the largest margin of 2.71%. This pattern is very similar to that under Condition 1, but land use is less altered in Condition 3, and the expansion trend is less than in Condition 1—suggesting that this pattern is likely to be subject to inefficient land use, as evident from the lowest scores for conversion cost for this simulation.

The experiment results demonstrate that our model provides decision makers and managers with spatially pleasing, efficient, and less conflicting land uses within a regional development setting, and that the proposed mutation modification significantly enhances overall fitness whilst also adhering to minimum patch size requirements. Furthermore, near-optimality is reached relatively quickly.

4. Discussion

In past decades, considerable progress has been made in multi-objective optimization using genetic algorithms for land-use optimization. However, most studies created impractical random initial patterns without minimum patch size constraint, resulting in the use of considerable computing toward a slow development of large patches engulfing small patches. In addition, social and spatial factors were only considered in objectives rather than in the implemented algorithm, resulting in differences between the optimal solution and actual land-use planning. Moreover, local optimum traps may appear because of inadequate internal process in GA causing premature convergence to non-optimal solutions.

To solve these problems, we developed the LIr-MSO model, which improves the initialization of patches, mutation, and elite strategy. First, patch size constraints for different land use was applied in initialization to obtain practical spatial agglomerations. Second, spatial allocation with implicit social relationships was enhanced through a two-step mutation operator to limit fragmentation and cell conversion. Finally, prematurity was minimized through a modified elite strategy.

The Dapeng Community was used as the study area to demonstrate the model’s ability to provide reliable support for land-use planning and sustainable development. The model used conversion cost, GDP, ESV, compactness, and conflict degree as objectives. According to the minimal block scale and relevant standards of urban planning in China, we set the minimum generated area of residential, commercial, industrial, and public land as 2, 0.5, 0.5, and 0.5 ha, respectively. A two-step mutation operator was used, which included block-based mutation and point-based mutation: block-based mutation used transfer coefficients to determine the change of land-use type with the lowest conversion cost for cells in a 3 × 3 window, while point-based mutation eliminated isolated cells by evaluating whether the center cell’s land-use type has the lowest count within a 3*3 window. In addition, elite strategy was improved by dividing the replacement mode into three conditions.

As the optimal land-use pattern is a tradeoff for all objectives, we found that Dapeng Community still has great potential for development. In the optimized pattern, residential, commercial, and industrial land increased by 7.81%, 107.07%, and 24.39% respectively. Residential land grew new space around the original area; commercial lands more than doubled through clustering and merging; arable land remained almost unchanged; and some commercial and residential land began to appear around industrial land in response to favorable alternative land-use transfer coefficients.

To demonstrate the priority of the model, experiments without improved operators were conducted with three experiments: without initialization or mutation or without both. The comprehensive fitness of the LIr-MSO model was 1.74%, 2.05%, and 2.71% higher than that of each condition, respectively. The experiments without the improved mutation had inefficient land use with less optimal conversion cost. This comparison shows that LIr-MSO can, under the assumed objectives and constraints, achieve an optimal more realistic land-use pattern. However, the land-use planning implications of the compactness objective being independent of patch size and in rounding off patches (due to rounding off reducing compactness) needs further consideration, especially in situations where existing transport networks cannot be redeveloped.

Lastly, the fitness trajectory showed a rapid initial phase to near optimality and a much slower parent–offspring mutation progress to the optimal pattern. Considering that the realistic simulation of millions of raster counts take hours or days to converge on a typical laptop, this observation could be used in realistic quick simulations that only run to near-optimal solutions to interactively assess the impacts of alternative nominalizations and weightings of objectives. This unexpected benefit of our model opens the possibility of interactive simulations involving planners and managers who can now assess the sensitivity of planning uncertainties before running a final convergent simulation.

5. Conclusions

Our proposed multi-objective genetic algorithm LIr-MSO model provides land-use planners with more realistic, efficient, and fitter land-use optimizations through much more detailed patterns, minimum initial patch size constraints, more balanced scaling, and novel genetic replacements in the mutation and elite strategy. Cell generation improvements result in a rapid increase of fitness to near-optimality, and the final fitness patterns appear to be practical with respect to the patch sizes of different land use and their placement relative to one another. The rapid fitness phase can facilitate quick realistic simulations that would allow planners to experiment with uncertainties and alternate models.

The reliability of the model was tested in an application to Dapeng Community, Shenzhen, which showed that the fitness of the model was 2.71% higher than that without improvements in the initialization and mutation. However, the percentage increase in fitness from initial to final pattern does not reflect substantial increases in some land use—for example, in commercial land use, which more than doubled.

However, improvements of the optimization model are still needed for practical applications. First, due to the complexity of urban development, more effective objectives should be considered. For example, the set of objectives should include increased pressures from population, traffic, accessibility, livability, and pollution control. The balance of these objectives in determining the final land use also needs further experimentation with planners to determine sets of suitable relative weights, which now seems computationally plausible using near-optimal solutions. Rapid near-optimal solutions would also allow decision makers to be intimately involved in land-use planning. Thus, it may be possible to combine models that can simulate the decision-making process from different agents—as in multi-agent systems—especially those involved in transport planning in view of potential tradeoffs between compactness and changes to transport networks. In addition, planning scenarios in the form of alternate regional development policies and plans should be designed to meet the different needs in practical land-use planning. For example, with the ecology-first concept, spatial optimization based on ecological security can be tested.