Reconstruction and Pattern Analysis of Historical Urbanization of Pre-Modern China in the 1910s Using Topographic Maps and the GIS-ESDA Model: A Case Study in Zhejiang Province, China

The pattern of urban land use and the level of urbanization in China’s pre-modernization period are of great significance for land use and land cover change (LUCC) research. The purpose of this study is to construct a 1910s spatial dataset of provincial land urbanization in pre-modern China. Using historical topographic maps, this study quantitatively reconstructs the built-up area of various cities in Zhejiang Province in the 1910s. The research indicates that: (1) During the early period of the Republic of China, there were a total of 252 cities and towns in Zhejiang Province, including 75 cities at or above the county level, 21 acropolis, and 156 towns. The total built-up area was 140.590 km2. (2) The county-level urbanization level had significant agglomeration characteristics. The overall urbanization rate of land was 0.135%. (3) Hot spots analysis showed that the Hang-Jia-Hu-Shao plain is hot spot. (4) The correlation coefficient between the city wall perimeter data recorded in the local chronicles and the measured city wall perimeter was 0.908. The research showed that the military topographic maps possessed a good application prospect for the reconstruction of urbanization levels. The research results provide direct evidence for urbanization and urban land use in China’s pre-modernization period.


Introduction
Land use and land cover change (LUCC) is one of the primary ways human activities affect the earth system [1], and it has complex regional and global ecological environmental impacts and responses. Among the various land cover changes, urbanized areas are the most severe form. Urbanization is the process of transforming agriculturally oriented areas into non-agricultural forms by changing the land use structures [2,3]. Unreasonable urbanization may cause serious ecological and environmental problems [4], such as the reduction of forests [5] and arable land [6], and the destruction of the rivers

Study Area
The study area of this article is Zhejiang Province (27°02′-31°11′ N, 118°01′-123°10′ E), which is located in the southeast coastal area of China. The total area is 105,500 km 2 and the terrain is sloping from southwest to northeast, with various landform types. Zhejiang Province is located in the middle of the subtropical zone and belongs to a monsoon humid climate with superior natural conditions produced during the Republic of China, an exploratory spatial data analysis (ESDA) in the ArcGIS platform is used to perform a spatial analysis and area calculation. The research results can provide historically basic data for the study of urbanization in pre-modern China. In addition, the urbanization data series can be extended by approximately a century.

Study Area
The study area of this article is Zhejiang Province (27°02′-31°11′ N, 118°01′-123°10′ E), which is located in the southeast coastal area of China. The total area is 105,500 km 2 and the terrain is sloping from southwest to northeast, with various landform types. Zhejiang Province is located in the middle of the subtropical zone and belongs to a monsoon humid climate with superior natural conditions

Study Area
The study area of this article is Zhejiang Province (27 • 02 -31 • 11 N, 118 • 01 -123 • 10 E), which is located in the southeast coastal area of China. The total area is 105,500 km 2 and the terrain is sloping Sustainability 2020, 12, 537 4 of 18 from southwest to northeast, with various landform types. Zhejiang Province is located in the middle of the subtropical zone and belongs to a monsoon humid climate with superior natural conditions [48]. The civilization of this area has been developed since ancient times and belongs to the traditional developed area in China [49].
In general, cities are large settlements formed by the aggregation of non-agricultural industries and non-agricultural population, including both small towns and large cities [50]. According to the Chinese traditional city classification criteria [18], it can be divided into a provincial capital city, prefecture-level cities, county-level cities, and other towns. Zhejiang Province is divided into 12 prefectures (Figure 3), and each prefecture is divided into several counties. The political center of Zhejiang Province is Hangzhou Fu, which belongs to the prefecture level, and the provincial capital city is Hangxian county. Other prefectures include Jiaxing, Huzhou, Jinhua, Quzhou, Yanzhou, Ningbo, Shaoxing, Wenzhou, Taizhou, Chuzhou, and Dinghai. There are 75 county-level administrative units, and the scope of jurisdiction is basically the same as that of modern Zhejiang Province. The boundaries of the administrative units at the prefecture and county levels were drawn based on the CHGIS-1911 dataset jointly developed by Fudan University and Harvard University [51]. Except for cities above the county level, the remaining urbanized areas are mainly towns, which provided commercial services to rural areas. It is worth noting that there is a special type of acropolis in the town of Zhejiang Province. This type of acropolis was originally a military castle, and gradually developed into a commercial town [52].
Sustainability 2020, 12, 537 4 of 20 [48]. The civilization of this area has been developed since ancient times and belongs to the traditional developed area in China [49]. In general, cities are large settlements formed by the aggregation of non-agricultural industries and non-agricultural population, including both small towns and large cities [50]. According to the Chinese traditional city classification criteria [18], it can be divided into a provincial capital city, prefecture-level cities, county-level cities, and other towns. Zhejiang Province is divided into 12 prefectures (Figure 3), and each prefecture is divided into several counties. The political center of Zhejiang Province is Hangzhou Fu, which belongs to the prefecture level, and the provincial capital city is Hangxian county. Other prefectures include Jiaxing, Huzhou, Jinhua, Quzhou, Yanzhou, Ningbo, Shaoxing, Wenzhou, Taizhou, Chuzhou, and Dinghai. There are 75 county-level administrative units, and the scope of jurisdiction is basically the same as that of modern Zhejiang Province. The boundaries of the administrative units at the prefecture and county levels were drawn based on the CHGIS-1911 dataset jointly developed by Fudan University and Harvard University [51]. Except for cities above the county level, the remaining urbanized areas are mainly towns, which provided commercial services to rural areas. It is worth noting that there is a special type of acropolis in the town of Zhejiang Province. This type of acropolis was originally a military castle, and gradually developed into a commercial town [52].

Military Topographic Map
The historical maps used in this study were primarily the military topographic maps of Zhejiang Province during the period of the Republic of China archived in the Key Laboratory of the Poyang Lake Wetland and Watershed Research Ministry of Education, Jiangxi Normal University. This set of topographic maps primarily consists of large-scale maps of 1:50,000. The mapping time was from the Republic of China year two (AD 1913) to the Republic of China year 25 (AD 1936), covering the entire province of Zhejiang and including parts of the Jiangxi, Anhui, and Jiangsu Provinces. The surveying and mapping organization was the Zhejiang Army Survey Bureau during the early stage. It was later changed to the Zhejiang Land Survey Bureau of the Staff Headquarters of the National Government. There are 345 maps of the whole province, and the size of each picture is approximately 46 cm × 35 cm. The longitude span of each map is 15′, and the latitude span is 10′. The elevation system assumes that the altitude of the Ziwei Garden in the old provincial Administration of Zhejiang is 50 m. The map schema adopts the original topographic map of the Republic of China year three (AD 1914). For a small portion of the missing area, topographic maps collected by the former National Government's

Military Topographic Map
The historical maps used in this study were primarily the military topographic maps of Zhejiang Province during the period of the Republic of China archived in the Key Laboratory of the Poyang Lake Wetland and Watershed Research Ministry of Education, Jiangxi Normal University. This set of topographic maps primarily consists of large-scale maps of 1:50,000. The mapping time was from the Republic of China year two (AD 1913) to the Republic of China year 25 (AD 1936), covering the entire province of Zhejiang and including parts of the Jiangxi, Anhui, and Jiangsu Provinces. The surveying and mapping organization was the Zhejiang Army Survey Bureau during the early stage. It was later changed to the Zhejiang Land Survey Bureau of the Staff Headquarters of the National Government. There are 345 maps of the whole province, and the size of each picture is approximately 46 cm × 35 cm. The longitude span of each map is 15 , and the latitude span is 10 . The elevation system assumes that the altitude of the Ziwei Garden in the old provincial Administration of Zhejiang is 50 m. The map schema adopts the original topographic map of the Republic of China year three (AD 1914). For a small Sustainability 2020, 12, 537 5 of 18 portion of the missing area, topographic maps collected by the former National Government's Ministry of the Interior organized by the Institute of the Modern History Academia Sinica were used as a supplement [53].
Another source of military topographic maps used in this study was a 1950s aerial survey map from the U.S. Army with a scale of 1:250,000. During the Anti-Japanese War, the Sino-US Cooperative Aerial Survey Team [54] conducted detailed aerial photogrammetry of eastern China. In the early 1950s, the Army Map Service of the US Corps of Engineers produced a set of 1:250,000 topographic maps of eastern China based on aerial survey data [55]. These are now archived in the University of Texas Library Pavilion [56]. This set of maps has detailed latitude and longitude information and projection properties that were used for spatial correction and information extraction during the early stage. Compared with this set of 1:50,000 military topographic maps of the Republic of China, aerial survey maps have better spatial coordinate information, so they can use the same place name points for spatial referencing and improve accuracy. Although the scales of the two sets of topographic maps are different, we extract all map information in the ArcGIS 10.2 platform based on the results of geographic referencing, and guarantee the reliability of the results through fine-tuning and compared with remote sensing imagery.

Spatiotemporal Scope
Since this study only involved the urban areas and towns in this set of topographic maps, the period statistics of the maps without relevant town information in the map were not counted. The mapping time of the related maps was primarily concentrated between 1913 and 1936 (Figure 4a), with a median of 1917 ( Figure 4b). Among them, the highest frequency occurred in 1916, accounting for approximately 25%. Therefore, it can be considered that this set of military topographic maps basically reflects the urban and town construction information of Zhejiang Province during the early period of the Republic of China. The scope of this study is Zhejiang Province during the early period of the Republic of China. Ministry of the Interior organized by the Institute of the Modern History Academia Sinica were used as a supplement [53]. Another source of military topographic maps used in this study was a 1950s aerial survey map from the U.S. Army with a scale of 1:250,000. During the Anti-Japanese War, the Sino-US Cooperative Aerial Survey Team [54] conducted detailed aerial photogrammetry of eastern China. In the early 1950s, the Army Map Service of the US Corps of Engineers produced a set of 1:250,000 topographic maps of eastern China based on aerial survey data [55]. These are now archived in the University of Texas Library Pavilion [56]. This set of maps has detailed latitude and longitude information and projection properties that were used for spatial correction and information extraction during the early stage. Compared with this set of 1:50,000 military topographic maps of the Republic of China, aerial survey maps have better spatial coordinate information, so they can use the same place name points for spatial referencing and improve accuracy. Although the scales of the two sets of topographic maps are different, we extract all map information in the ArcGIS 10.2 platform based on the results of geographic referencing, and guarantee the reliability of the results through fine-tuning and compared with remote sensing imagery.

Spatiotemporal Scope
Since this study only involved the urban areas and towns in this set of topographic maps, the period statistics of the maps without relevant town information in the map were not counted. The mapping time of the related maps was primarily concentrated between 1913 and 1936 (Figure 4a), with a median of 1917 ( Figure 4b). Among them, the highest frequency occurred in 1916, accounting for approximately 25%. Therefore, it can be considered that this set of military topographic maps basically reflects the urban and town construction information of Zhejiang Province during the early period of the Republic of China. The scope of this study is Zhejiang Province during the early period of the Republic of China.

Process Flow and the GIS-ESDA Model
The research process of this study primarily included the following three steps ( Figure 5). First, based on the detailed spatial coordinate information of the 1950s aerial survey map from the U.S. Army, a preliminary georeferencing of the set of 1:50,000 military topographic maps from the 1910s was performed. Then the position of the preliminary located topographic map was fine-tuned based on the remote sensing image to ensure the accuracy of the spatial information. Finally, based on the Zipf rule and the GIS-ESDA model, a spatial analysis and area measurement of the urban information on the historical topographic map was completed.

Process Flow and the GIS-ESDA Model
The research process of this study primarily included the following three steps ( Figure 5). First, based on the detailed spatial coordinate information of the 1950s aerial survey map from the U.S. Army, a preliminary georeferencing of the set of 1:50,000 military topographic maps from the 1910s was performed. Then the position of the preliminary located topographic map was fine-tuned based on the remote sensing image to ensure the accuracy of the spatial information. Finally, based on the Zipf rule and the GIS-ESDA model, a spatial analysis and area measurement of the urban information on the historical topographic map was completed.

Preliminary Georeferencing
The spatial positioning of this set of topographic maps primarily included two steps. First, this set of topographic maps were initially registered and positioned in the ArcGIS 10.2 platform according to the stitching relationship. Second, based on the space coordinate information of the aerial map from the US military, the military map of the 1910s was spatially located according to some same place names and mountains or rivers.

Fine-Tuning of Topographic Map
Fine-tuning of topographic map primarily included two steps. First, modern Landsat satellite TM/ETM+ remote sensing images [57] were used to repeatedly fine-tune the position of the topographic map based on obvious natural landmarks (such as mountains and rivers). Second, control points were utilized for overall error control and evaluation.

Spatial Analysis Zipf Distribution Rule
It is generally believed that the size distribution of cities in the different levels in a region follows a certain law. G. K. Zipf proposed the law in 1949: where Si is the scale of the i-th city; S0 is the scale of the first city; and Ri is the rank of the i-th city. In view of the fact that the Zipf distribution law is an ideal state of order-scale distribution in cities, it is generally expressed in terms of the Rotka model [58,59]:

Preliminary Georeferencing
The spatial positioning of this set of topographic maps primarily included two steps. First, this set of topographic maps were initially registered and positioned in the ArcGIS 10.2 platform according to the stitching relationship. Second, based on the space coordinate information of the aerial map from the US military, the military map of the 1910s was spatially located according to some same place names and mountains or rivers.

Fine-Tuning of Topographic Map
Fine-tuning of topographic map primarily included two steps. First, modern Landsat satellite TM/ETM+ remote sensing images [57] were used to repeatedly fine-tune the position of the topographic map based on obvious natural landmarks (such as mountains and rivers). Second, control points were utilized for overall error control and evaluation.

Spatial Analysis Zipf Distribution Rule
It is generally believed that the size distribution of cities in the different levels in a region follows a certain law. G. K. Zipf proposed the law in 1949: where S i is the scale of the i-th city; S 0 is the scale of the first city; and R i is the rank of the i-th city.
In view of the fact that the Zipf distribution law is an ideal state of order-scale distribution in cities, it is generally expressed in terms of the Rotka model [58,59]: where q is the fractal dimension of the urban system. Generally, the logarithm of both sides of Equation (2) can be taken to obtain: Then the fractal dimension, q, is obtained using the method of least squares. Studies [60] have shown that when q > 1, the distribution of urban systems in the region is more concentrated, and large cities are prominent. When q < 1, the distribution of urban systems in the region is relatively scattered, and small or medium cities develop.
Exploratory Spatial Data Analysis, the ESDA The ESDA is a collection of spatial data analysis technologies. It primarily reveals spatial correlation characteristics and patterns by studying the spatial heterogeneity of data and the visualization of spatial distribution patterns [61]. The ESDA mainly includes a global spatial autocorrelation, a local spatial autocorrelation, a hotspot analysis, and a spatial trend analysis.
This study used the Moran's I index for the global spatial autocorrelation analysis [62,63]. The larger the Moran's I index, the higher the degree of spatial autocorrelation in the region. The calculation method is: where n is the number of sample points; y i and y j are the attribute values of the sample points; and w ij is the spatial weight matrix. The local spatial autocorrelation was identified using the Getis-Ord G i * index [62,63]. The G i * index is normalized, and the Z value is obtained. A higher Z value indicates that this area is a hot spot gathering area, while a lower Z value indicates that this area is a cold spot gathering area. The calculation method for the Z value of the Getis-Ord G i * index is: where n is the number of sample points; x i and y j are the attribute values of the sample points; and w ij is the spatial weight matrix.

Kernel Density Analysis
Kernel density analysis is a non-parametric test method that can be used to analyze the density of spatial point distributions. The basic principle is to estimate the theoretical distribution of sample points in a region using a kernel density function and to convert the density of discrete sample points into density values that are continuously distributed in space. A kernel density analysis can identify the concentrated area of a spatial point feature distribution; that is, the hotspot distribution area. The calculation method is: where k (·) is the kernel function; r is the analysis radius; and x − x i is the distance between the point x to be estimated and the sample point x i .

GeoDa and ArcGIS Platforms
The calculation process involved in this study was mainly performed using GeoDa 0.9.5 and ArcGIS 10.2 software platforms. The main calculation indicators include the Moran's I Index, the Getis-Ord G i * index, and the kernel density value.

Reconstruction of Urban Built up Areas
Based on the measured topographic maps of the urban areas, the spatial distribution ( Figure 6a) and built-up area (Figure 6b) of each city or town were obtained using the ArcGIS 10.2 platform. Statistics show that during the period of the Republic of China (Table 1), there were a total of 252 cities and towns of various types in the Zhejiang Province, including 75 cities at the county level and above, 21 acropolises, and 156 towns. The total built-up area of towns and cities in the province was 140.590 km 2 , of which the total built-up area of cities at the county level and above was 91.764 km 2 , with an average value of 1.223 km 2 . The total built-up area of acropolises was 8.051 km 2 , with an average value of 0.383 km 2 . The total built-up area of towns was 40.775 km 2 , with an average value of 0.261 km 2 . Of all of the cities and towns in Zhejiang Province, the largest built-up area was Hangxian, which is the capital city of the province, with an area of 14.778 km 2 . The smallest was Maqiao town, Haining County, with an area of only 0.023 km 2 .
Sustainability 2020, 12, 537 8 of 20 GeoDa and ArcGIS Platforms The calculation process involved in this study was mainly performed using GeoDa 0.9.5 and ArcGIS 10.2 software platforms. The main calculation indicators include the Moran's I Index, the Getis-Ord Gi * index, and the kernel density value.

Reconstruction of Urban Built up Areas
Based on the measured topographic maps of the urban areas, the spatial distribution ( Figure 6a) and built-up area (Figure 6b) of each city or town were obtained using the ArcGIS 10.2 platform. Statistics show that during the period of the Republic of China (Table 1)

Spatial Pattern of Urbanization
The level of urbanization can be measured from the perspectives of urbanization of land and urbanization of population. Because urban population data for historical periods are difficult to obtain, this study used the ratio of the total area of urban built-up areas and land area at the county level scale as urbanization rate.

Spatial Pattern of Urbanization
The level of urbanization can be measured from the perspectives of urbanization of land and urbanization of population. Because urban population data for historical periods are difficult to obtain, this study used the ratio of the total area of urban built-up areas and land area at the county level scale as urbanization rate.
With the support of GeoDa 0.9.5 and ArcGIS 10.2, the Moran's I index Equation (7) of the county-level urbanization in Zhejiang Province during the 1910s was 0.296 (Z-Score = 4.992, p < 0.001). This indicated that its spatial distribution had significant agglomeration characteristics (Figure 7). More precisely, the high-value area of the urbanization level was adjacent to the high-value area, and the low-value area was adjacent to the low-value area.
With the support of GeoDa 0.9.5 and ArcGIS 10.2, the Moran's I index Equation (7) of the countylevel urbanization in Zhejiang Province during the 1910s was 0.296 (Z-Score = 4.992, p < 0.001). This indicated that its spatial distribution had significant agglomeration characteristics (Figure 7). More precisely, the high-value area of the urbanization level was adjacent to the high-value area, and the low-value area was adjacent to the low-value area. As can be seen from Figure 8a, among all county-level units, Hangxian County had the highest urbanization rate at 1.463%. Jingning County had the lowest at 0.014%, and the overall urbanization rate of Zhejiang Province was 0.135%. Further analysis of the local spatial autocorrelation cold and hot spots (Figure 8b) showed that the Z-score of the Getis-Ord Gi * index Equation (5) had obvious high-high/low-low aggregation characteristics in the spatial distribution. More precisely, the hot spots were primarily distributed in the Hang-Jia-Hu-Shao plain area nearby Hangzhou Bay, and the cold spot concentration areas were near west Zhejiang and southwest Zhejiang. The area near Wenzhou Bay belonged to the medium distribution area of the Z-score value. The spatial distribution of the cold and hot spots were generally in the direction of northeast-southwest; that is northeast Zhejiang Province belonged to hotspots and the southwest belonged to cold spots. As can be seen from Figure 8a, among all county-level units, Hangxian County had the highest urbanization rate at 1.463%. Jingning County had the lowest at 0.014%, and the overall urbanization rate of Zhejiang Province was 0.135%. Further analysis of the local spatial autocorrelation cold and hot spots (Figure 8b) showed that the Z-score of the Getis-Ord G i * index Equation (5) had obvious high-high/low-low aggregation characteristics in the spatial distribution. More precisely, the hot spots were primarily distributed in the Hang-Jia-Hu-Shao plain area nearby Hangzhou Bay, and the cold spot concentration areas were near west Zhejiang and southwest Zhejiang. The area near Wenzhou Bay belonged to the medium distribution area of the Z-score value. The spatial distribution of the cold and hot spots were generally in the direction of northeast-southwest; that is northeast Zhejiang Province belonged to hotspots and the southwest belonged to cold spots.

Urban System Structure
The expressions of the rank-scale distribution of all of the 252 cities and towns in Zhejiang Province during the Republic of China were obtained using the Zipf distribution rule (Equations (1)-

Urban System Structure
The expressions of the rank-scale distribution of all of the 252 cities and towns in Zhejiang Province during the Republic of China were obtained using the Zipf distribution rule (Equations (1)-(3)): y = −1.010x + 3.352 (r 2 = 0.903, p < 0.001). According to Equation (3), we can see that the fractal dimension is q = 1.010 >1. It can be seen that, as a whole, the distribution of the Zhejiang urban system was relatively concentrated during the 1910s. Furthermore, the kernel density analysis method Equation (6) in ArcGIS 10.2 was used to obtain the spatial pattern of all cities and towns in Zhejiang Province, such as cities at the county level and above, towns, and acropolises. From Figure 9a, it can be seen that the kernel density values of spatial distribution of all of the 252 cities and towns in Zhejiang Province during 1910s formed two high-value central areas on the north and south sides of Hangzhou Bay. In addition, two secondary nuclear density regions existed near the Ning-shao Plain and Wenzhou Bay. From Figure 9b, it can be clearly seen that the spatial distribution of 75 county-level cities and above in Zhejiang Province is evenly distributed, and no relatively high-value areas had been formed throughout the province. Since cities at the county level and above were basically political institutions in pre-modern China, their even spatial distribution is conducive to administrative control.  Nevertheless, the spatial density of the county-level cities and above in northern Zhejiang was still higher than that in southern Zhejiang, which is consistent with the fact that the northern counties in Zhejiang are more distributed than the southern ones. According to the spatial distribution kernel density maps of 156 towns in the province (Figure 9c), it can be seen that the spatial distribution pattern is basically consistent with the spatial structure of all cities and towns in Zhejiang Province (Figure 9a). Moreover, a high value range of kernel density is primarily distributed on the both sides of Hangzhou Bay and areas such as the Ning-Shao Plain. This also shows that the urban system pattern of Zhejiang Province during the Republic of China was primarily determined by the spatial distribution of towns. The areas, such as the Hang-Jia-Hu Plain and the Ning-Shao Plain, have always developed market economies, so a large number of towns have developed under county-level cities. Due to historical reasons, such as resistance to pirating and coastal defenses, Zhejiang's coastal areas Nevertheless, the spatial density of the county-level cities and above in northern Zhejiang was still higher than that in southern Zhejiang, which is consistent with the fact that the northern counties in Zhejiang are more distributed than the southern ones. According to the spatial distribution kernel density maps of 156 towns in the province (Figure 9c), it can be seen that the spatial distribution pattern is basically consistent with the spatial structure of all cities and towns in Zhejiang Province (Figure 9a). Moreover, a high value range of kernel density is primarily distributed on the both sides of Hangzhou Bay and areas such as the Ning-Shao Plain. This also shows that the urban system pattern of Zhejiang Province during the Republic of China was primarily determined by the spatial distribution of towns. The areas, such as the Hang-Jia-Hu Plain and the Ning-Shao Plain, have always developed market economies, so a large number of towns have developed under county-level cities. Due to historical reasons, such as resistance to pirating and coastal defenses, Zhejiang's coastal areas also have had a large number of acropolises. From Figure 9d, it can be seen that the areas with high kernel density values in the urban spatial distribution of the acropolis are primarily concentrated along the coast of eastern Zhejiang and the southern shore of Hangzhou Bay. In addition, another high kernel density area formed in the Wenzhou Bay area.

Trend Analysis
From an analysis of the spatial trend of the built-up areas of all of the 252 cities and towns (Figure 10a), it can be seen that the distribution trend rises from west to east and rises from south to north. This also shows that the spatial distribution of all the urban areas in Zhejiang Province has a trend of being high in the northeast and low in the southwest. The spatial trend analysis of the urbanization level of the 75 county-level regions in Zhejiang Province also indicates this result (Figure 10b).
Sustainability 2020, 12,537 12 of 20 This also shows that the spatial distribution of all the urban areas in Zhejiang Province has a trend of being high in the northeast and low in the southwest. The spatial trend analysis of the urbanization level of the 75 county-level regions in Zhejiang Province also indicates this result (Figure 10b).

Uncertainty Analysis
Although the military topographic map of the Republic of China used in the reconstruction of

Uncertainty Analysis
Although the military topographic map of the Republic of China used in the reconstruction of urban built-up areas in this study was completed using modern surveying and mapping technology, the technological conditions and cartographic levels at that time had historical limitations. Therefore, it is necessary to evaluate the errors of these topographic maps and perform an uncertainty analysis. By referring to the principle of the "point-to-point" test [26,42], 26 relevant control points were randomly selected in the military topographic map after secondary correction and fine-tuning according to the principle of uniform distribution throughout the province. These were then compared with the corresponding position in the TM/ETM+ remote sensing image. The error range between the control points on the remote sensing image and the points on the topographic map were calculated in ArcGIS 10.2 (Table 2), and error range distribution map was interpolated by kriging method in ArcGIS 10.2 ( Figure 11). The error range of the 26 control points was between 0.142 and 1.348 km, and the average error was 0.553 km. Considering that these sets of military topographic maps had relatively small errors after secondary correction and projection fine-tuning, they can be used for the reconstruction of relevant geographical elements in this area.

Comparison of the Perimeter of the City Wall and Historical Records
During the historical period, various national archives and chronicles have generally been recorded in the shape of the city and the perimeter of the city wall, but they are often described non-

Comparison of the Perimeter of the City Wall and Historical Records
During the historical period, various national archives and chronicles have generally been recorded in the shape of the city and the perimeter of the city wall, but they are often described non-quantitatively in the form of text or hand-drawn maps. For example, during the Qing Dynasty (AD 1644-1911) in the Emperor Jiaqing Restoration Unification Book, the city wall form of Fuyang County in Zhejiang Province was recorded as "The perimeter is six li. Four gate. There is a ring moat around the city. Southeast of the city near the river . . . " He et al. [64] made full use of the city wall perimeter data recorded in the historical books and obtained the city area by considering the shape of the city as a square or a circle. They also obtained a set of built-up areas of cities in the 18 provinces of the Qing Dynasty. Since the construction of the city wall was basically not done in Zhejiang Province after the period of the Republic of China, the city information reflected on the topographic map of the Republic of China can be considered to represent the perimeter. Therefore, it can be considered the actual city wall shape and area of the city during the Qing Dynasty.
Based on the city wall perimeter data recorded in the local chronicles and in the Emperor Jiaqing Restoration Unification Book, the city wall perimeter information of 55 cities above the county level was obtained, including 12 prefecture-level cities. In addition, the unit of measurement of the length of the Qing Dynasty li was converted according to 1 li = 0.576 km [65], and the area formula is: where A i is the area of the i-th city; and C i is the circumference of the i-th city. The city wall perimeter was converted into the area (km 2 ). By comparing the historical city wall perimeter data and city wall perimeter data based on the topographic map of the Republic of China (Figure 12a), it was found that the correlation coefficient between the two groups reached 0.908 (p < 0.001). It can be seen that, although the city perimeter data recorded in the historical books are limited in terms of the approximate number and low precision, the overall length is significantly consistent with the actual city perimeter length reconstructed from the military topographic map. The results of the regression analysis showed that the city perimeter data recorded in the historical books was approximately 20% higher than the measured value. It should be noted that, since Hangxian is the capital of Zhejiang Province, it has a large wall perimeter and area than other cities. After excluding Hangxian data from the model, the correlation coefficient is 0.834 (p < 0.001).
although the city perimeter data recorded in the historical books are limited in terms of the approximate number and low precision, the overall length is significantly consistent with the actual city perimeter length reconstructed from the military topographic map. The results of the regression analysis showed that the city perimeter data recorded in the historical books was approximately 20% higher than the measured value. It should be noted that, since Hangxian is the capital of Zhejiang Province, it has a large wall perimeter and area than other cities. After excluding Hangxian data from the model, the correlation coefficient is 0.834 (p < 0.001). Similarly, there is a significant correlation between the measured city area and the city area calculated using Equation (7) (Figure 12b), and the correlation coefficient is 0.955 (p < 0.001). The area calculated using Equation (7) is approximately 80% higher than the measured area. The correlation coefficient is 0.843 (p < 0.001) after excluding Hangxian data from the model. Similarly, there is a significant correlation between the measured city area and the city area calculated using Equation (7) (Figure 12b), and the correlation coefficient is 0.955 (p < 0.001). The area calculated using Equation (7) is approximately 80% higher than the measured area. The correlation coefficient is 0.843 (p < 0.001) after excluding Hangxian data from the model.

Comparison of Land Urbanization and Population Urbanization
Land urbanization and population urbanization are two aspects of the urbanization process. As detailed statistics of urban land and population in history are often subject to many limitations, the relationship between the two groups needs to be further examined. A regression analysis of the urban area of the county and the total population of the county in Zhejiang Province for the middle period of the Republic of China showed that the correlation coefficient between them was 0.792 (p < 0.001) (Figure 13a). It can be seen that there was a significant correlation between the urban area of the county and the total population. Generally, the area of the county with a large population has a larger urban built-up area. After excluding Hangxian data from the model, the correlation coefficient is 0.810 (p < 0.001). The reason for this change may be that Hangxian is the capital city of Zhejiang Province and has a high proportion of the urban population in the total population.

Comparison of Land Urbanization and Population Urbanization
Land urbanization and population urbanization are two aspects of the urbanization process. As detailed statistics of urban land and population in history are often subject to many limitations, the relationship between the two groups needs to be further examined. A regression analysis of the urban area of the county and the total population of the county in Zhejiang Province for the middle period of the Republic of China showed that the correlation coefficient between them was 0.792 (p < 0.001) (Figure 13a). It can be seen that there was a significant correlation between the urban area of the county and the total population. Generally, the area of the county with a large population has a larger urban built-up area. After excluding Hangxian data from the model, the correlation coefficient is 0.810 (p < 0.001). The reason for this change may be that Hangxian is the capital city of Zhejiang Province and has a high proportion of the urban population in the total population. Although there was no direct urban demographic data in the Internal Affairs Survey Statistics Form, the county town population and general town population can be used to approximate the urban population in a certain county. In addition, the population urbanization level of the county can be obtained using this method. The calculation results showed that the overall population urbanization rate in Zhejiang Province was 13.49%, with the highest value of Hangxian being 64.71% and the Although there was no direct urban demographic data in the Internal Affairs Survey Statistics Form, the county town population and general town population can be used to approximate the urban population in a certain county. In addition, the population urbanization level of the county can be obtained using this method. The calculation results showed that the overall population urbanization rate in Zhejiang Province was 13.49%, with the highest value of Hangxian being 64.71% and the lowest value of Qingyuan being 1.08%. The regression analysis of the land urbanization level and the population urbanization level showed that the correlation coefficient between them was 0.799 (p < 0.001) (Figure 13b). It should be noted that the reason for the high urbanization rate of the population in Hangxian County might be because the original data included the entire population of the provincial capital, but it did not distinguish between urban and agricultural populations. After excluding Hangxian data from the model, the correlation coefficient is 0.609 (p < 0.001). Despite this fact, the coupling relationship between land urbanization and population urbanization can still be studied using the two sets of data.

Conclusions and Prospects
(1) During the period of the Republic of China, there were a total of 252 cities and towns in Zhejiang Province including 75 cities at or above the county level, 21 cities in acropolis, and 156 towns. The total built-up area of towns and cities in the province was 140.590 km 2 , of which the total built-up area of cities at the county level and above was 91.764 km 2 , with an average value of 1.223 km 2 . The total built-up area of the acropolis was 8.051 km 2 , with an average value of 0.383 km 2 . In addition, the total built-up area of towns was 40.775 km 2 , with an average value of 0.261 km 2 .
(2) During the period of the Republic of China, the county-level urbanization level of Zhejiang Province had significant agglomeration characteristics. Among the county-level units, the urbanization rate was highest in Hangxian County at 1.463%, and lowest in Jingning County at 0.014%. The province's overall land urbanization rate was 0.135%. An analysis of the cold and hot spots showed that the hot spots were primarily distributed in the Hang-Jia-Hu-Shao area near Hangzhou Bay. The cold spots were primarily concentrated in the western and southwestern Zhejiang. The region near Wenzhou Bay belonged to the medium distribution area of the Z value.
(3) According to the Zipf distribution rule, it was concluded that the distribution of Zhejiang's urban system during the Republic of China belonged to the centralized distribution type. The results of the kernel density analysis showed that the urban spatial distribution formed two high-value central areas on the north and south sides of Hangzhou Bay and two secondary high-density areas near the Ning-Shao Plain and Wenzhou Bay.
(4) The correlation coefficient between the city wall perimeter data recorded in the local chronicles and the measured city wall perimeter was 0.908, showing that the historical books and chronicles had high reliability. However, the results of the regression analysis showed that the city perimeter data recorded in the historical books was approximately 20% higher than the measured value.
(5) The province's overall urbanization rate was 13.49%, of which the highest value of Hangxian County was 64.71% and the lowest value of Qingyuan County was 1.08%. During the period of the Republic of China, there was a significant correlation between the urbanization rate of county land and the urbanization rate of the population, and the correlation coefficient was 0.799.
This study systematically demonstrated the relationship between the city wall perimeter data in historical local chronicles and the measured urban built-up area using the 1:50,000 large-scale topographic map data during the Republic of China, which was in pre-modern China. Additionally, according to the regression analysis, the relational equation of the historical city wall perimeter and the measured urban area was obtained. Most parts of China have relatively systematic historical records of the city wall perimeter that can be converted to the urban area using the relationship formula, thus providing basic data for a long-term urbanization level sequence study. In addition, by performing a comparison with the population data of Zhejiang Province during the Republic of China, the correlation between the land urbanization and the population urbanization data was obtained. This correlation provides an empirical research result for understanding the relationship between different measures of urbanization during historical periods.
However, due to the accuracy limitations of historical data, the urban built-up area reconstructed in this paper can only represent the situation in the early and middle period of the Republic of China. The coupling relationship between land urbanization and population urbanization is only a preliminary discussion. Later land-use reconstruction work requires in-depth research on the spatial projection correction and the control of reconstruction errors to ensure the reliability of reconstruction data. Additionally, the spatial and temporal resolution of land use reconstruction for an historical period should be further improved, and the dynamic pattern evolution analysis of land use and urbanization should be performed on a more detailed scale.