An Automatic Generalization Method for a Dense Road Network Area Considering Spatial Structural Features as Constraints

Abstract

Road networks are the skeletal elements of topographic maps at different scales, and road selection is a prerequisite for implementing continuous multiscale spatial representations of road networks. The mesh-based approach is a common, advanced and powerful method for road selection in dense road areas in which small meshes are removed and road segments with the least importance in each mesh are eliminated. However, small meshes in a map can be classified into two types: aggregated small meshes and isolated small meshes. The number of the former is small, and that of the latter is large. Existing methods are generally applicable for the latter, and some or even most spatial characteristics will be lost when they are applied to the former; as a result, the road selection quality will be affected. Therefore, as a supplement to the mesh-based selection method, this paper proposed an automatic generalization method of dense road network areas (areas formed by aggregated small meshes) considering spatial structural features as constraints. First, the aggregated areas of small meshes were identified based on the number of adjoining small meshes, and the boundaries of aggregated areas are extracted and used as hard constraints during mesh elimination. Second, the starting meshes were redefined by simultaneously considering the edge features and mesh density of small meshes, and an ordinal elimination algorithm was proposed to eliminate the meshes in the stroke connection direction. Third, road selection was implemented by identifying the starting meshes and sequentially processing the related mesh pairs. This iterative process continued until all mesh densities of the newly formed meshes are beyond the threshold or the problem becomes a simple elimination problem involving two adjoining small meshes or one isolated small mesh. Finally, a 1:10,000 standard topographic road map for Jiangsu Province, China, was used for validation. The experimental results showed that in the aggregated areas with two small meshes, 31% of the areas obtained the same selection results by using the mesh-based method and the proposed method, and the remaining 69% obtained a more compact result with the proposed method. Moreover, for all aggregated areas with more than two small meshes, the spatial distribution structure of small meshes was preserved better by the proposed method.

1. Introduction

Road networks are the skeletal elements of topographic maps at different scales. To obtain continuous multiscale spatial representations, road networks require generalization [1,2,3]. The key and core step in the generalization processing is road selection, which aims to reduce the level of detail in a road network by retaining important roads and eliminating less important roads while also maintaining the main road spatial characteristics and network structure [4,5,6], as well as essential semantic, topological, and geometrical properties [7]. With these selection constraints, maintaining the spatial distribution structure of the road network is the most difficult part; hence, the structured selection of road networks has been recently studied in detail. Among these studies, the mesh-based approach [8] has been reported to be a common, advanced, and powerful model to complete road selection work in dense road areas [5,9,10]. The research in this paper also falls into this field.

To obtain a satisfactory selection result for road networks, an enormous amount of research has been carried out worldwide. Existing methods can be divided into three categories, graph-theory-based methods, stroke-based methods, and mesh-based methods. For graph-theory-based methods, the topologic structure of a road network is captured by a graph structure [11,12]. However, the graph structures are not conducive to maintaining the connectivity and integrity of roads [13]. Thomson et al. [14,15] proposed a stroke-based method by considering the good continuation principle in Gestalt psychology. Road segments in a stroke are treated as a unit in the selection process; that is, they are simultaneously deleted or selected. By considering the importance of each road, the stroke-based method performs selective omission by eliminating the least-important road stroke iteratively [16,17]. The connectivity and integrity of roads are preserved well in this method, and it was widely used by scholars [18,19]. However, the geometric patterns of a road network, which are usually represented as mesh patterns of road segments, cannot be reflected in this method. To address the limitations of stroke-based methods, Chen et al. [8] proposed a road selection method based on the mesh density of a road network, which is the aforementioned mesh-based method. In this method, all small meshes with mesh densities beyond a threshold value are removed in sequence, and road segments with the least importance in each small mesh are eliminated. This method is able to effectively reflect the overall network characteristics and local density variations across different regions. Li and Zhou [20] showed that the mesh-based approach performs best for areal pattern roads (i.e., road segments that span one or two meshes). Wu et al. [21] and Touya [5] also applied the mesh-based approach to road selection in urban areas, obtaining satisfactory results.

Based on the number of adjoining small meshes, the road network can be divided into three types: (1) isolated small meshes, which have no adjoining small meshes, as illustrated in the gray part in Figure 1; (2) aggregated small meshes, which have two or more adjoining small meshes, as illustrated in the blue part in Figure 1; and (3) other meshes, which are not small meshes, as illustrated in the yellow part in Figure 1. Existing studies are more practical and applicable to areas with isolated small meshes, which account for the vast majority of all small meshes on a map. However, there are some aggregated small meshes in residential blocks or large factory areas. The typical spatial characteristics of these areas, which refer to the cluster boundary morphology, mesh arrangement and mesh uniformity, will be destroyed by using the existing mesh-based selection methods. Although few areas of aggregated small meshes generally exist on a map, these meshes still play an important role in the quality of the selection process. Therefore, as a supplement to the mesh-based selection method, an automatic generalization method of dense road network areas considering spatial structural features as constraints was proposed in this study.

The paper is organized as follows. Section 2 presents the advanced mesh-based road network selection method and analyzes its limitations for road network generalization in dense road network areas. Section 3 introduces an automatic generalization method for dense road network areas considering spatial structural features as constraints. Section 4 provides a series of experiments to validate the reliability and superiority of the proposed method. The conclusions are presented in Section 5.

2. Related Works

2.1. Current Mesh-Based Selection Methods

The traditional mesh-based approach was proposed by Chen et al. [8], in which mesh is defined as a closed region surrounded by several road segments. The detailed steps in this method are given as follows [8,10]:

Step 1: Construct the node-arc-polygon topology of the road network and identify meshes based on topological polygons;

Step 2: Calculate the density of each mesh in a road network according to Equation (1):

$$D=\frac{P}{A}$$

(1)

where D is the mesh density, P is the perimeter of the mesh, and A is the area of the mesh;

Step 3: Determine the mesh density threshold according to the target scale by using a statistical algorithm based on sample data;

Step 4: Identify the meshes that have a density beyond the threshold, define them as small meshes, and add them to the candidate set to be processed;

Step 5: Sort the small meshes in the candidate set in descending order and handle the one with the largest mesh density first;

Step 6: Compare the importance of bounding road segments and eliminate the least-important segment. The importance of each edge segment is calculated with Equation (2) [8,10]:

$$I=\max \left\{{\omega }_1{C,\omega }_2{D}_S{,\omega }_3{L}_S{,\omega }_{4}L\right\}$$

(2)

where I indicates the importance index of each segment, C is the normalized road class, D_s is the normalized degree of the road stroke where the segments are located, L_s is the normalized length of the road stroke, L is the normalized length of the road segment, max{*} is a function to determine the importance of the segments, and ${\omega }_i$. (i = 1, 2, 3, 4) represents the weight of each parameter. Before the calculation, these four parameters are normalized to (0, 1) through the min-max normalization algorithm, and ${\omega }_1$, ${\omega }_2$, ${\omega }_3$, and ${\omega }_4$ are set to 10⁴, 10³, 10², and 10, respectively;

Step 7: Merge the remaining segments with the adjoining mesh and form a new mesh;

Step 8: Calculate the density of the new mesh and add it to the candidate set if its mesh density exceeds the threshold;

Step 9: Repeat steps 5–8 until the set of small meshes is empty.

The above process is illustrated in Figure 2. In Figure 2a, the mesh with the largest density of 0.42 is first merged with the mesh with a density of 0.40, and the relatively unimportant segment L is eliminated. The density (0.16) of the new mesh is then updated (Figure 2b). The process is iteratively repeated until the densities of all the meshes are smaller than a given threshold, final elimination result is illustrated in Figure 2c.

The flowchart is shown in Figure 3.

2.2. Limitations of the Mesh-Based Selection Method

Small meshes in a map can be classified into two types: aggregated small meshes and isolated small meshes. Existing studies are generally practical and applicable for areas with small isolated meshes, which account for the vast majority of small meshes on a map. The selection results in these areas can preserve the spatial structure and density characteristics of road networks well by existing methods. However, there are a small number of areas with aggregated small meshes in residential blocks or large factories region, and the typical spatial characteristics of these areas will be destroyed by using the traditional mesh-based selection method. Although the number of aggregated small meshes on a map is small, the selection results of these areas still affect the quality of the whole area.

As shown in Figure 4a, the aggregated area includes 12 small meshes, which are marked by black numbers. Figure 4b displays the selection result of the mesh-based method, and the red numbers describe the elimination order. The regular grid structure formed by these 12 small meshes is not captured through this elimination order, and the spatial structure characteristics of the original data are obviously changed. In contrast, the human cartographer with a global perspective seeks to obtain an output in which the main characteristics of these small meshes are abstracted well [4]. As shown in Figure 4c, this result provides a better understanding of the road network structure.

3. Methodology

As a supplement to the common mesh-based selection method, an automatic generalization method of dense road network areas considering spatial structural features as constraints is proposed in this paper. This method consists of the following four steps: (1) identifying aggregated areas for small meshes and the extracting their boundaries to establish hard constraints during road selection; (2) redefining starting meshes for mesh elimination by simultaneously considering the boundary features and mesh density; (3) establishing an ordinal elimination algorithm, in which meshes are eliminated in the stroke connection direction; and (4) iteratively processing all small meshes inside the aggregated areas until the mesh densities of the newly formed meshes are beyond the threshold or the problem is transformed into the simple elimination problem involving two adjoining small meshes or one isolated small mesh.

3.1. Identification of Aggregated Areas for Small Meshes and the Extraction of Their Boundaries

Determination of the mesh density threshold is a core step in identifying small meshes. Chen et al. [8] proposed a sample statistical algorithm in which the appropriate threshold value is determined based on the split nodes of the mesh density distribution curves at two scales before and after the road network generalization. For example, practical experiments proved that 0.016 m/m² is considered an appropriate value of the density threshold for transformation from a 1:10,000 to 1:50,000 urban road network data [10]. Correspondingly, meshes with mesh densities larger than the threshold are defined as small meshes. This paper focuses on road selection operation in aggregated areas of small meshes. The adjoining relationship between two small meshes is calculated based on the node-arc-polygon topology. If two small meshes (i.e., polygons) share a common segment, they are considered adjoining to each other.

Previous studies have shown that the boundary feature is an important parameter that can describe the spatial structure of an object in maps with large scale. In this paper, the boundary of each aggregated area of small meshes is extracted according to the topology and used as a hard constraint, i.e., the boundary is directly preserved and is not eliminated during the elimination process.

3.2. Redefining Starting Meshes for Mesh Elimination

For the common mesh-based method, the small mesh with the largest mesh density is selected as the starting mesh. However, if the starting mesh is located in the middle of an aggregated area of small meshes, the mesh elimination result is difficult to control, and it is easy to change the spatial distribution structure of the original road network. Hence, we redefine the starting mesh by simultaneously considering the boundary features and mesh density of small meshes, the starting mesh is changed from one single small mesh to a pair of adjoining small meshes. The reason for this is as follows. First, compared with processing only one single mesh at a time, processing pairs of adjoining meshes is more conducive to controlling the local spatial structure of the road network. Second, road selection from the boundary areas to the inside of a network can avoid the problems caused by selection from the center of the road network. Considering the above two constraints, the starting meshes for mesh elimination are redefined in this paper as two adjoining meshes with the largest mesh density and located in the outermost part of the aggregated areas. If two or more mesh pairs have equal metrics, identify the common road segments between the adjoining small meshes in each mesh pair, and compare the importance of the road segments. The mesh pair with the least-importance road segment is taken as the starting mesh pair.

The sum of the mesh density of two adjoining meshes is easy to calculate, and the location of the meshes is determined by their first-order adjacent field. For a small mesh, its first-order adjacent field is the area formed by the small meshes that share common segments or nodes with it. Obviously, the small meshes located at the boundary have the lowest number of small meshes in the first-order adjacent field.

As depicted in Figure 5a, for the mesh with a density of 0.42 (lavender), the meshes in the first-order adjacent fields are those with densities of 0.40, 0.38 and 0.28 (pink), and the number of small meshes in the first-order adjacent field is 3. Similarly, in Figure 5b, for the mesh with a density of 0.38 (lavender), the first-order adjacent field includes the meshes with densities of 0.42, 0.40, 0.28, 0.26 and 0.34 (pink), and the number is 5. Hence, the mesh with a density of 0.42 is located closer to the network boundary than that with a density of 0.38. Moreover, the small mesh pair formed by the meshes with densities of 0.42 and 0.40 (orange) has the largest summed mesh density, and these meshes are selected as the starting meshes.

3.3. Ordinal Elimination Algorithm

As described in Section 2.2, the order of mesh elimination is another important factor for preserving the spatial structure of a road network. To obtain an improved elimination result, we propose an ordinal elimination algorithm in which the meshes are eliminated according to the stroke connection directions, i.e., road segments connection directions. In this study, strokes are constructed based on the arc-node topology of the road network with consideration of the names and directions of road segments at a junction. This process was described in detail in our previous study [22].

By using the ordinal elimination algorithm for a mesh pair selected to be processed in aggregated areas of small meshes, the common segment between the meshes is marked firstly. Then, the stroke where the common segment is located is identified, and the other mesh pairs that share other road segments in this stroke are captured. Unlike the common mesh-based method, this method uses an elimination process that is not in descending order based on the mesh density but is instead the stroke connection order. The advantage of this approach is that the loss of spatial structure caused by some legacy strokes can be avoided. As illustrated in Figure 6a, the meshes with densities of 0.42 and 0.40 (orange) are defined as starting meshes, and the stroke S (blue) where their common segment is located is marked (Figure 6b). The other mesh pairs that share other road segments in this stroke are eliminated in sequence, and the elimination result is shown in Figure 6c. Figure 6d illustrates the elimination result based on considering the mesh density in descending order; it is clear that the elimination result based on the stroke connection order can better abstract the original spatial structure.

3.4. Iterative Processing

To eliminate all aggregated areas of small meshes, an iterative processing algorithm is proposed. The workflow for this algorithm involves the following steps:

Step 1: Determine the mesh density threshold value and classify small meshes by using a statistical algorithm based on sample data;

Step 2: Identify the aggregated areas of small meshes from the road network and add them to the candidate set to be processed;

Step 3: Randomly select an aggregated area, extract the boundary, and use the result as a hard constraint during the mesh elimination process;

Step 4: Identify the starting meshes in this aggregated area based on the algorithm in Section 3.2 and handle the starting meshes first;

Step 5: Identify the common segment of starting meshes and detect the stroke associated with the common segment;

Step 6: Determine the spatial relationship between the stroke and the aggregated area. If the stroke crosses the aggregated area, splitting the aggregated area into two sub-aggregated areas based on this stroke and processing step 7. Otherwise, go to step 7 directly.

As shown in Figure 7a, road stroke S1 crosses the aggregated area, and road segments a and b are inside the aggregated area; if these two road segments are removed directly, the connectivity of S1 will be lost. Therefore, the following constraint should be used in the elimination process: if a road stroke crosses an aggregated area, split the aggregated area into two sub-aggregated areas based on this stroke, and implement the iterative processing algorithm for each sub-aggregated area. As shown in Figure 7b, the original aggregated area of small meshes is divided into two sub-aggregated areas (green and blue) based on stroke S1;

Step 7: Implement the ordinal elimination algorithm in Section 3.3, and finish processing the mesh pairs that connect the stroke;

Step 8: Identify the new starting meshes without considering the newly formed meshes obtained in steps 5–7.

As depicted in Figure 8a, the two adjoining small meshes in orange are the starting meshes in the first round, and the newly formed meshes (orange in Figure 8b) will not participate in the identification of the starting meshes in the second round. Therefore, as illustrated in Figure 8c, the two adjoining small meshes shown in blue are the starting meshes in the second round, and the final selection result is shown in Figure 8d;

Step 9: Repeat steps 5–8 until no new starting meshes can be identified; at this point, the first round of mesh elimination is finished;

Step 10: Calculate the density of all newly formed meshes and compare these values with the threshold to identify new small meshes in this aggregated area;

Step 11: Repeat steps 4–10 until the mesh densities of newly formed meshes all exceed the threshold or only two small meshes exist in this aggregated area.

Step 12: Eliminate the common segment of the two small meshes, calculate the density of the new mesh, and compare this value with the threshold; if the threshold is exceeded, proceed to the next step;

Step 13: Process this new small mesh by using our stroke edge-based method [10];

Step 14: Repeat steps 3–13 until all the aggregated areas in the candidate set have been processed.

The flowchart is shown in Figure 9.

4. Experiments and Analyses

4.1. Experimental Data and Environment

The proposed method in this paper is embedded into the WJ-III mapping workstation developed by the Chinese Academy of Surveying and Mapping. The experimental data come from a 1:10,000 standard topographic road map of Jiangsu Province in China. The spatial scope of the experimental data is 23.91 × 18.67 km², and the generalization target scale is 1:50,000. The mesh density threshold is set to 0.016 m/m² at a scale of 1:50,000. The original 1:10,000 experimental data are shown in Figure 10a. To verify the reliability and effectiveness of the proposed method, the traditional mesh-based method was used for comparative analysis, and the selection result of the traditional mesh-based method and our method is demonstrated in Figure 10b,c, respectively. The operating environment used to run the software system is an Intel Core I7-3770 CPU with the Windows 7 64-bit operating system, a main frequency of 3.2 GHz, 16 GB of memory, and a 1024 GB solid-state hard disk.

As illustrated in Figure 10, the selection result of the proposed method satisfied the specification for maps at 1:50,000 in China according to the China Standardization Office (GB/T 12343.1—2008) [23], in which it is required that the important and high-grade roads should be retained firstly, and attention should be paid to maintaining the density difference and shape characteristics of the road network.

4.2. Basic Statistical Analysis of the Experimental Data

There are 1782 meshes in the experimental data, including 471 small meshes and 1311 non-small meshes. The detailed classification of the 471 small meshes is shown in

Table 1. The detailed classification of small meshes.

Isolated Small Meshes	Aggregated Area of Small Meshes
112	64 aggregated areas (include 359 small meshes)
	Number of small meshes in aggregated area
	Min	Max	Average
	2	47	5.6

.

Table 1 shows that 23.77% of small meshes are isolated small meshes, and the remaining small meshes form 64 aggregated areas. Thus, the small meshes in the experimental area display obvious aggregation characteristics. In the aggregated areas of small meshes, the minimum number of small meshes is 2, which indicates that two adjoining small meshes are needed to form an aggregated area, and the largest aggregated area consists of 47 small meshes.

4.3. Reliability and Effectiveness Analysis

(1)

Comparative analysis of aggregated areas with two small meshes

For the 32 aggregated areas with two small meshes, the two different methods yielded 10 areas with the same selection result, accounting for 31.25% of all aggregated areas in this case; the least-important segment of each mesh is just the common segment, which is why the selection results were the same. The number of areas with selection results that were considerably different based on the two different methods was 22, accounting for 68.75% of all aggregated areas in this case. The reason for this result is that the least-important segment of each mesh is not always the common segment for these meshes. Shape similarity (SS) is used to quantitatively evaluate the selection results. The formula is as follows:

$$\mathit{SS}=\frac{{\sum }_{i=1}^n{\mathit{BArea}}_i}{{\sum }_{i=1}^n{\mathit{AArea}}_i}$$

(3)

where n represents the total number of meshes related to an aggregated area, i is the i^th small mesh, $BAre{a}_i$ indicates the area of the i^th small mesh before road selection, and $AAre{a}_i$ indicates the area of the mesh in which the original i^th small mesh is contained after road selection. Obviously, when the value of SS is 1, the shape similarity of road meshes before and after road selection is optimal.

The SS for each method in each aggregated area is calculated. To compare the results of the two methods and understand the degree of differences between them, the calculation results are shown in Figure 11.

As depicted in Figure 11, for the mesh-based method, the minimum value of SS is 0.01, which means that one of (or both of) the two small meshes in this area is merged into a very large adjoining mesh, and the area of the original small meshes changed dramatically. In contrast, the corresponding SS values of the two small meshes are 1 based on the proposed method; this result suggests that the two small meshes can form a new mesh with a mesh density that satisfies the threshold. A similar situation occurs for the other 21 aggregated areas. For the aggregated areas No. 23–No. 32, the SS values of the two methods are all 1, indicating that these two methods obtained the same selection results, and the shape features of these 10 areas are preserved well. As a consequence, the proposed method has a smaller influence on the spatial structure of the road network than does the traditional method, and the local features are better preserved. Figure 12 illustrates two representative areas: aggregated areas with the same selection results and aggregated areas with different results.

From Figure 12a,b, the selection results of the two different methods for the first typical area are the same, and the least-important segment of small mesh is just the common segment a. However, as depicted in Figure 12c,d, the least-important segment of the mesh is not always the common segment a; consequently, the two methods yield different selection results for the second typical area. Compared with the mesh-based method, the proposed method produces a selection result that better maintains the compact structure of small meshes and preserves the connectivity of the road stroke where road segment b is located.

(2)

Comparative analysis of aggregated areas with more than two small meshes

For the 32 aggregated areas with more than two small meshes, the shape similarity (SS) and a new index named the coefficient of variation of the area (CVA) are used to quantitatively evaluate the selection results. The proposal of CVA is based on the following considerations. Normally, the geometric area of each mesh in the aggregated areas is homogenous. When the scale changes from large to small, a reasonable generalization result should maintain this characteristic, which means that the geometric area of each mesh should not vary considerably.

The formula for the CVA is as follows:

$$\mathit{CVA}=\frac{{\sigma }_{\mathit{area}}}{{\mu }_{\mathit{area}}}$$

(4)

where ${\sigma }_{area}$ represents the standard deviation of the area of small meshes in the aggregated area and ${\mu }_{area}$ is the average area of small meshes in the aggregated area. The CVA is used to measure the uniformity of the geometric area of a small mesh in an aggregated area. The smaller the CVA value is, the more uniform the area distribution.

The SS and CVA of each method in each aggregated area are calculated, and the results are shown in Figure 13.

As depicted in Figure 13a, the SS values of the 32 aggregated areas obtained by the mesh-based method have a range of [0.03, 0.78], which indicates that the area of all small meshes change to varying degrees. In contrast, the SS values of these aggregated areas obtained by the proposed method are all 1, and boundary constraints ensure that the spatial structure of the original aggregation area does not change in general. Additionally, compared with the CVA values of the original aggregated areas, the values of 27 aggregated areas obtained by the proposed method decrease. In contrast, the values of 26 aggregated areas obtained by the mesh-based method increase, as depicted in Figure 13b. When the boundaries of aggregated areas are eliminated by the mesh-based method, the small meshes are merged to adjacent other meshes, which form larger meshes. This will increase the area difference between meshes, and the value of CVA will also increase. The comparison of CVA value curves shows that the density distribution of the selection results obtained by the proposed method is more uniform than that of the results of the mesh-based method.

For the 32 aggregated areas with more than two small meshes, the selection results of the two methods are different. The number of small meshes included in these 32 aggregated areas lies within the range of (3, 47). Hence, this paper selects the aggregated areas with the fewest (3), median number (8) and maximum number (47) of small meshes as typical regions to visually compare the generalization results of the two different methods shown in Figure 14.

As seen in Figure 14, the differences in the three typical regions for aggregated areas with more than two small meshes are as follows.

(1)

For the aggregated area with the minimum number of small meshes, the mesh-based method eliminates road segments a, b and c based on their lower importance levels, leading to the change of the spatial structure of the road network in this area. In contrast, the proposed method eliminates only isolated road strokes a and d by considering the boundary constraint, and a more reasonable selection result that better reflects the aim of a cartographer;

(2)

For the aggregated area with the median number of small meshes, the mesh-based method eliminates road segments a, b, c, d, h and g, which are located in four different strokes, the spatial structure in this area is inevitably changed. In contrast, the proposed method eliminates only two road strokes (one is formed by segments a and b, and the other is formed by segments e, f and g) inside the aggregated area using the ordinal elimination algorithm. Consequently, the selection results maintain the original distribution structure well on the whole;

(3)

For the aggregated area with the maximum number of small meshes, the method proposed in this paper better maintains the mesh distribution density than that based on the traditional method. The obvious difference in the selection results of these two methods occurs in the area where road strokes S1, S2 and S3 are located. It is clear that the preservation of stroke S1 has an important influence on the uniformity of the distribution density. Similarly, the preservation of stroke S3 yields a clear outline. Moreover, stroke S2 is completely preserved by the proposed method, which reflects the selection intention of human cartographers.

5. Conclusions

The traditional mesh-based selection method is a common, advanced, and powerful road network generalization method in which the small meshes are sequentially removed and the least-important road segment in each mesh is eliminated. However, for areas with aggregated small meshes that encompass residential blocks or large factories, the corresponding typical spatial characteristics are easily destroyed. Therefore, as a supplement to the mesh-based selection method, this paper proposes an automatic generalization method for dense road network areas considering spatial structural features as constraints. The following conclusions were drawn from the experimental validation and analysis using actual data.

(1)

The experimental data include 1782 meshes, approximately 74% of the meshes are non-small meshes that do not need to be processed during road selection; correspondingly, the remaining 26% of meshes are small meshes. Approximately 24% of the small meshes are isolated small meshes that can be reasonably eliminated by the existing method;

(2)

Approximately 14% of the small meshes form aggregated areas with 2 small meshes, and the two different methods yield 31% areas with the same selection results. The selection results of the remaining areas are considerably different, and the proposed method obtain a compact aggregate that better reflects the local spatial structure of the road network compared to the results of the traditional method;

(3)

Approximately 62% of the small meshes form aggregated areas with more than two small meshes, the proposed method performs better in maintaining the density difference and shape characteristics of the road network. Through the comparison of the CVA curves, the mesh density distribution of 90% of the selection results obtained by the proposed method is more uniform than the original data.

The refined starting meshes and the new elimination order in our method avoid the random error that occurs with the traditional method and are conducive to maintaining the original road network structure. However, our method has three limitations. First, the map scale in which this method can be applied is limited, this method is more suitable for the road network selection in large scale (>1:100,000). In a small scale, the detailed spatial structure feature is not necessary to be preserved, and the aggregated areas may be eliminated as a whole. Second, the spatial structure considered in our paper comes from a global perspective, and finer spatial structures cannot be captured and distinguished. Third, the selection threshold is a fixed value in the experimental area. If this value can be adaptively calculated according to the spatial characteristics of the local area, the selection result will be more reasonable. Hence, future research should focus on two topics. First, different aggregation patterns, such as grid patterns, narrow patterns, and radial patterns, should be systematically summarized and classified because the fine division of aggregation patterns is conducive to maintaining the spatial distribution characteristics of road networks. Second, the mesh density threshold plays an important role in road selection, and the setting of an appropriate threshold in road selection for multiscale and multitype road network datasets requires further study.