Reconstructing and tracing the evolution of the road networks in the Haidai region of China during the Bronze and Early Iron Ages

Abstract

Reconstructing ancient transportation networks is critical to studying past human mobility patterns. China’s Haidai region was a thriving political and economic hub during the Bronze and Early Iron Ages. We used GIS spatial analysis techniques to build a “Settlement Interaction Model” based on archaeological data from the Haidai region during the Bronze and Early Iron Age (Shang Dynasty, Western Zhou Dynasty, Spring & Autumn Period, and Warring States period). The eight-level road network maps with traffic attributes were distinguished based on topography and settlement size. The total lengths of the road networks were estimated to be 19,112 km in the Shang Dynasty, 35,269 km in the Western Zhou Dynasty, 51,555 km in the Spring & Autumn Period, and 77,456 km in the Warring States Period, with the average road flows of 6.6, 31.7, 42.8, and 75.5, respectively. The Z score and one-sample t-test (p < 0.01) confirmed the reliability of the reconstructed road networks. The Shang Dynasty saw the sporadic appearance of simple road routes. More complex routes emerged during the Western Zhou Dynasty and Spring & Autumn Period. The road networks were finally built during the Warring States Period. The development of road networks was closely related to population growth and urbanization. Exploring methods for reconstructing road networks may help us uncover ancient road networks and better understand ancient cultural exchanges.

Introduction

Cultural exchange and integration facilitated social progress and civilization development, and the communication network, particularly the road network, gradually formed and evolved as a result, promoting human migration, material exchange, and cultural communication. Primitive paths were formed alongside human dispersal during the Paleolithic Age, including the African path of Homo sapiens migration [1,2,3]. The establishment of initial exchange routes then had an impact on the spread of language and culture in the Neolithic-Bronze Age [4,5,6,7], especially with food globalization and the spread of metal smelting technology [8,9,10,11]. Increased food production and productivity have significantly aided population growth and social development. The “Silk Road” stretched across Eurasia for nearly two millennia [12, 13], and the “New Voyages” opened in the fifteenth-sixteenth centuries, connecting the entire world and contributing to social development and civilizational collisions [14].

Simulating ancient road networks is crucial for understanding ancient communication [15,16,17,18]. Archaeologists and historians used historical documents and archaeological data to infer the discontinuous transmission routes of materials and technologies in ancient societies, such as crop diffusion [19], domesticated animals [20] and metal smelting [11, 21]. Nowadays, with the aid of powerful spatial analysis tools like Geographic Information Systems (GIS), continuous communication routes or networks have been simulated, such as the early human diffusion routes and modes [1, 22, 23], mobile network of animal husbandry groups [12], ancient maritime exchange route network [24, 25], and the ancient routes from historical records (e.g., the Silk Road and the Tangbo ancient route) [26, 27]. Additionally, road flow simulation was used to better understand the characteristics of ancient road networks. Howey [28] took the road network as a circuit, where the current intensity in the circuit was depicted by the probability of passing through a node. In addition, a “From Everywhere to Everywhere” (FETE) model was used to simulate the natural formation process of roads without requiring input from origins and destinations [29]. A high-resolution mobility network in Asian highland grasslands was created by simulating pastoralists’ seasonal nomadic activities [12]. According to these research ideas, spatial analysis techniques should be used to reconstruct road networks at various levels using known archaeological data in an area with frequent human activity, as well as to recover complex communication patterns and road network evolution patterns among complex settlements within cultural zones to gain a better understanding of past inter-population communication.

The Haidai region was an important birthplace of Chinese civilization and one of the political and economic centers from the Bronze to the Early Iron Age [30, 31]. During the Bronze to Early Iron Age, cities and populations increased rapidly [32, 33], and the Haidai region advanced significantly before being incorporated into the Shang kingdom (1600–1046 B.C.). This region saw the emergence of regional centers such as Daxinzhuang, Subutun, and Qianzhangda, as well as the discovery of numerous bronze artifacts and oracle bone inscriptions [34, 35]. The “System of Enfeoffment” then triggered the peak of city building in the Haidai region during the early Zhou Dynasty [36], with Qi, Lu, and Ju states forming and their capitals gradually becoming regional metropolises. For example, Qufu was one of the first cities and the longest-serving capital [37], while Linzi was one of the largest at the time, covering an area of more than 15 km² [38]. Furthermore, commercial activities in the Haidai region, such as salt production and iron smelting, boosted human communication. Linzi, as the capital of Qi State, was a major iron smelting center during the Warring States Period [39], and a large number of salt workshop site clusters were also investigated in the Xiaoqing River in north Shandong [40]. The complex trade network gradually developed, aided by tremendous advancements in craftsmanship and business organization patterns in the Haidai region during the Bronze and Early Iron Ages [41,42,43].

The Settlement Interaction Model (SIM) was developed to trace the traffic road network in the Haidai region, and it is based on the concept of least-cost paths (LCPs) and urban interaction theory in economic geography [44]. Taking into account settlement size and topography, river systems, vegetation, and other natural factors, road networks were graded according to the intensity of interaction caused by road flow between settlements. The road networks in the Haidai region were then reconstructed for four periods, namely the Shang Dynasty (1600–1046 B.C.), Western Zhou Dynasty (1046–771 B.C.), Spring & Autumn Period (771–476 B.C.), and Warring States Period (476–221 B.C.), with lengths and morphologies of the road networks calculated and analyzed in each period. Furthermore, the influencing factors were briefly discussed. This approach provided a new perspective on regional mobility, which is critical in revealing the evolutionary process of road networks and investigating the relationship between human communication and the physical environment.

Study area

The Haidai region (114°48′–122°42′E, 34°23′–38°17′N) is located in the east part of North China, which covers an area between the Bohai Sea (Hai) and Mount Tai (Dai). The Haidai cultural area gradually formed around 4000 BC and belonged to the Neolithic Age [30, 31]. It is centered on the Taiyi Mountain system, which stretches east to the Yellow Sea and west to the Cangzhou-Dezhou-Puyang line, with the Bohai Sea to the north and the North Jiangsu Plain to the south. The Haidai region is characterized by mountainous terrain with plains and basins interspersed. The complex topography includes mountains in the center, low-lying areas in the southwest and northwest, and gently undulating hills in the east. The Taiyi Mountain Range, Lushan Mountain, Mengshan Mountain, and other mountains make up the majority of Shandong’s mountainous area. Mount Tai has an elevation of 1500 m above sea level at its peak, with a relative height of 1300 m between the top and the foot of the mountain. As a result, the rugged mountainous landscape of central Shandong hampered north–south communication. Furthermore, the Yellow River flows from southwest to northeast into the Bohai Sea, and the river valleys formed by rivers in the Haidai region cut through the mountain systems, likely aiding human migration through the Haidai region.

Furthermore, the Haidai region is located in the warm temperate monsoon climate zone, with an average annual temperature of 11–14 °C and an average annual precipitation of 550–950 mm. Shandong’s vegetation is primarily made up of coniferous forests and deciduous broad-leaved forests, which reflect the region’s complex and diverse topography and the characteristics of a warm temperate monsoon climate.

Materials

Many archaeological sites’ data were primarily derived from “The Atlas of Chinese Cultural Relics”, which included volumes on Shandong, Jiangsu, Anhui, Hebei, and Henan [45,46,47,48,49], as well as published excavation reports. Each site’s details were meticulously recorded, including its number, era, area, latitude, and longitude. The digital elevation model (DEM) was obtained from a joint publication of the Ministry of Economy, Trade, and Industry of Japan (METI) and the United States National Aeronautics and Space Administration (NASA), specifically the 30 m SRTM (Shuttle Imaging Radar Topographic Mapping) (http://www.gscloud.cn), and a DEM was obtained using the “mosaic” and “mask extraction” functions in ArcGIS Pro. Furthermore, the ground cover data was based on a 1:4,000,000 vegetation map (http://westdc.westgis.ac.cn), while the original data for the water systems came from the National Geographic Information Resource Catalogue Service (www.webmap.cn). To achieve more accurate ancient topography, ancient maps were reregistered using volume 1 of “The Historical Atlas of China” and related elements, such as ancient lakes and rivers, were digitized in the Haidai region [50] (Fig. 1).

Fig. 1

Locations of the Haidai region of China

Methods

In spatial archaeology, computer simulations have become valuable tools for reconstructing ancient paths and road networks [17, 51]. The LCPs calculation method is widely used in archaeology, particularly for reconstructing and analyzing ancient pathways and migrated networks [1, 15, 17]. However, these methods face numerous challenges. The technique is primarily helpful for calculating a single optimal path, which presents a limitation when investigating multiple potential routes between settlements. Most LCPs methods prioritize cost considerations in modeling, ignoring the importance of specific settlements and the intensity of their connections. As a result, networks reconstructed via simple LCPs approaches often lack hierarchical structure and centrality, while exhibit high redundancy [52]. Reconstructing intricate communication networks for a large number of settlements (for instance, exceeding 1000) poses a significant challenge, especially when attempting to reveal the complexity of interactions across a broader geographical scale. To address these issues, we used the SIM in our research, a novel approach combining LCPs, settlement interactions, and the comparison of traffic to hydrologic flow. Figure 2 shows the overall flowchart for this method. Initially, the Pandolf equation was used to calculate the friction surface for movement, referred to as the travel metabolic rate. The accumulated cost surfaces were then calculated using the friction surface and settlement points. Incorporating past road influences into the accumulated cost surfaces is useful for studying the evolution of road networks [53]. However, using incomplete archaeological records as data input can complicate road simulations and lead to biased estimates. Past roads can only be included as a parameter in SIM if sufficient archaeological data is obtained in the future. This study proposes that the significance of roads within the network, referred to as “traffic flow” in later sections, is determined by the interaction strength between connected settlements, which is primarily influenced by distance and size [54]. Given the scarcity of archaeological data, a parameter-free radiation model was used to estimate the intensity of ancient exchanges [55]. The final roads were classified as surface drainage and extracted using a surface hydrological analysis method. Subsequently, a method for assessing the accuracy of the road network was implemented. Finally, the distinction between “road” and “route” should be clarified to avoid misunderstandings. A “road” was defined as a pathway with specific traffic flow characteristics, the magnitude of which was determined by the intensity of interactions between settlements along the way. A “route” was defined as a set of road segments that connect one settlement to another, emphasizing the continuity of space and direction.

Fig. 2

Settlement Interaction Model workflow diagram

Data collection and pre-processing

Extensive data collection was carried out before restoring the road network to ensure the accuracy of road network generation. Ancient archaeological sites were collected between this region and the surrounding areas, including Shandong and relevant areas in Hebei, Henan, Anhui, and Jiangsu provinces, to rebuild the entire road network. The text descriptions of sites were gathered from the atlas and reports, as described in Sect. “Materials”. Then, all of these locations were manually marked on Google Earth and their decimal coordinate information was entered into the database.

The “Slope” tool in ArcGIS Pro calculates slope data from the DEM, which serves as the terrain foundation for model operation and analysis. However, it is easier for route simulation to traverse lakes and rivers because its lower slopes are typically depicted as flat terrain in lakes and rivers based on the DEM, resulting in low friction in route simulation and guiding simulated routes across these lake-stream systems. Previous research assigned a preset friction value to the lake-stream system to block the route and avoid “puddle jumping” [15, 56]. Fábrega-Álvarez and Parcero-Oubiña [57] found that assigning a slope value of 15° to all lake-stream systems improves route modeling accuracy.

In addition, land cover should be considered, as it significantly impacts walking energy consumption [15, 29]. Due to a lack of high-resolution paleovegetation data in this area, modern vegetation cover data was found to be an acceptable compromise [12, 29, 56]. As a result, this study relied on the People’s Republic of China Vegetation Map for surface cover data. To account for varying resistance effects, assign a reasonable terrain factor [η in Eq. (1)] to each vegetation type [29, 58]. The forest factor was set to 1.5, shrub and farmland to 1.2, and water system to 1.8. The topographic factors were then assigned to different vegetation grids using ArcGIS Pro’s “Raster Calculator.”

Modeling of LCPs friction surface based on Pandolf equation

LCPs essentially represent the impact of geographic constraints on travel within the spatial landscape as “resistance”—the “friction surface” proposed in this paper. This is commonly measured in terms of expenditure on time or energy [29, 59,60,61]. In this study, the Pandolf equation was used to establish a more quantifiable measure of resistance within the spatial landscape, which translated some geographic constraints into travel costs, such as terrain, slope, and walker fundamentals. This equation was used to calculate the unit energy cost of movement in a study area because it is helpful for calculating energy expenditure while slow walking [60]. The standard equation read as follows:

$$MRT=1.5\cdot w+2\cdot {\left(\frac{l}{w}\right)}^{2}+\eta \cdot (w+l)\cdot \left(1.5\cdot {v}^{2}+0.35\cdot {v}_{{\text{e}}}\cdot {s}_{{\text{p}}}\right)$$

(1)

where MRT is the metabolic rate of travel (in W or J/s), w is the mass of the traveler (in kg), l is the weight of the load (in kg), η is the terrain factor, v is the walking speed (in m/s), and s_p is the percentage slope. The average height of men was around 165 cm in the Yellow River basin of China during the Shang and Zhou periods [62], and the average weight of ancient people was set at 60 kg (w), in consideration of the nutritional status of humans at that time, which should be poorer than the modern people. The load weight was set at 5 kg (l) because White and Barber [29] demonstrated that the load weight has less effect on the main routes, and a smaller load is better for identifying high-traffic routes. In addition, the walking speed v was derived from the well-known Tobler hiking function [61], which portrayed the nonlinear relationship between human walking speed and slope, expressed as:

$$v=\left(6\cdot {e}^{-3.5\cdot \left|k+0.05\right|}\right)\cdot \left(\frac{1}{3.6}\right)$$

(2)

$$k = {\text{tan}}\left( {\frac{{\theta \pi }}{{180}}} \right)$$

(3)

where v denotes the walking speed (in m/s), k denotes the percentage slope, and θ is the slope size (in degrees). Finally, a conversion factor was multiplied to convert the original unit of velocity from km/h to m/s. To obtain the percentage slope S_p, a conversion of θ was required:

$${s}_{{\text{p}}}=\frac{\theta }{90^\circ }\cdot 100$$

(4)

After calculating the MRT magnitude, it was multiplied by a unit conversion factor (1W = 0.000239 kcal/s) to convert it into resistance (in kcal/s), and divided it by the previously calculated velocity (v) to obtain the resistance value per unit distance with the following equation:

$$kcal=\frac{(0.000239)MRT}{v}$$

(5)

where kcal represents the energy required to travel a unit distance (in kcal/m). This value was entered into the ArcGIS Pro “Distance Accumulation” tool to calculate the accumulated cost surfaces, representing the energy consumption from the origin to any point. Consequently, this supported the computation of the route between two settlements.

Evaluating settlement interaction intensity and weighted communication distance

The settlement interaction intensity is the scale of population exchanges, which is the number of people traveling along the roads that connect these settlements. The Radiation Model (TRM) was adopted as a quantitative indicator to quantify the intensity of interaction [55]. This formula is consistent with the first law of geography and an important archaeology hypothesis: the closer two settlements are geographic, the more frequently they interact [44, 63,64,65]. Notably, this model does not rely on observational data for parameter calibration. The specific formulation is as follows:

$${T}_{ij}\equiv {T}_{i}{p}_{ij}={T}_{i}\frac{{m}_{i}{n}_{j}}{\left({m}_{i}+{s}_{ij}\right)\left({m}_{i}+{n}_{j}+{s}_{ij}\right)}$$

(6)

T_ij represents the average flux from location i to location j, while T_i denotes the number of travelers departing from location i. The probability is denoted by p_ij from location i to location j, and m_i and n_j represent the population of location i and j respectively. Additionally, s_ij is the total population within a circle of radius r_ij centered on location i, excluding m_i and n_j. Hence, the interaction magnitude between two settlements can be described accurately by inputting the population information and settlement locations. The resulting T_ij serves as a measure of the strength of the inter-cluster interaction and was utilized as a flow weighting factor of the road network.

The official population count in China began during the Han Dynasty (202 B.C.–220 A.D.), so prior population data can only be derived from historical estimates provided by scholars. Previous research found that Shandong Province had populations of approximately 1.1 million, 1.4 million, 2 million, and 4.2 million during the Shang, Western Zhou, Spring & Autumn, and Warring States periods, respectively [66]. Furthermore, archaeological site areas served as proxies for ancient populations, and conversion coefficients were used to convert settlement areas into population numbers [67, 68]. As a result, settlement areas in the Haidai region served as proxies for populations during the Bronze and Early Iron Ages.

Beyond the scope of the SIM inclusion section, the method outlined below aims to better understand the distance variations associated with inter-settlement communication, a weighted average of the settlements’ external communication distance was calculated based on the intensity of settlement interaction:

$${c}_{i}=\frac{\sum_{j=1}^{k}{T}_{ij}{r}_{ij}}{\sum_{j=1}^{k}{T}_{ij}}$$

(7)

Here, c_i (in km) represents the radius index of cluster exchange, k indicates the total number of settlements other than i in a particular period, and T_ij and r_ij have the same meaning as above. In General, larger weighted averages of communication distance suggest a preference for long-distance communication, whereas smaller averages suggest a preference for short-distance communication.

Extraction of road network based on the flow accumulation model

In addition to communication routes, traffic flow should reflect the importance of roads. When ancient human mobility is envisioned as water flowing across the surface in spatial terms, roads exhibit a continuous convergence from high-cost to low-cost areas, analogous to tributaries merging and forming a larger mainstream [12]. The flow accumulation process in surface analysis was used to achieve the aforementioned simulative process. As a result, a hydrologic flow accumulation algorithm was used within the SIM on the accumulated cost surfaces to express road semantics comprehensively [69]. In particular, the accumulated cost surfaces correspond to the surface raster in the hydrologic accumulation algorithm, which significantly impacts the potential directions of human mobility at each location. Every cell in the space is assigned a flow weight factor, determined by the interaction intensity between the cell and the origin settlement. As a result, weights designated as T_ij are assigned to cells at other settlement locations using Eq. (6), whereas cells that do not correspond to settlement locations have zero weights (Fig. 3). As a result, the distributions of rivers represent the morphology of ancient roads, while river flow reflects the number of people who use them. This algorithm was applied to each iterative process of the site point set, and the resulting flow rasters were aggregated to generate the road network. It should also be noted that traffic flow is only considered a dimensionless quantity due to the unknown number of people [70]. It represents a proportional relationship and can be interpreted as a relative number of people rather than an absolute number.

Fig. 3

Schematic diagram of the flow accumulation method (flow rasters overlaid to form a grid road network) for Sites A, B and C

Road network accuracy evaluation

Our evaluation method is based on a similar assessment approach to that of Frachetti et al. [12], and it validates the simulation’s accuracy by quantifying the spatial correspondence between independent historical site points (hereinafter referred to as control points) and SIM’s high-flow roads. However, good spatial correspondence does not imply that the results are non-random. Drawing on the work of Frachetti et al. [12], statistical analysis was used to compare the spatial correspondence between control points or random points and the road network in order to determine the significance of their differences.

Initially, control points adjacent to key routes within the study area are identified, independent of settlement points used in road simulations, with an assumed number of control points of n. Then, a 500 m buffer zone was established around these control points to evaluate the adjacency relationship between routes and control points within the buffer zone. The maximum traffic flow is then extracted within the buffer zone, representing the road with the highest flow adjacent to the control point. A total flow value is calculated by adding all extracted flow values.

The null hypothesis assumes no significant difference exists between the proportion of control points and random points on high-flow roads. To simulate spatial correspondence due to randomness, n random points (equal to the number of control points) are placed within the study area, and the total flow value is calculated for each round of random iteration using the same steps. This random iteration process is repeated 200 times, yielding 200 random point samples. The distribution of total flow values from random points is then compared to the count of total flow values recorded from control points, resulting in the Z Score of the single-sample T-test. If the simulated network is sufficiently accurate, the Z Score will have an absolute value greater than 2, indicating a significant deviation from the random model. Simultaneously, the P-value from the single-sample T-test should be less than 0.05, indicating a statistically significant difference and allowing us to reject the null hypothesis.

Results

The “Settlement’s interaction model” was executed by the model builder tool in ArcGIS Pro. The raster-represented traffic road networks were created by iterating over each point in each period and adding them together to make the final ones over four periods. We use the Natural Breaks method and manual adjustments to ensure statistical significance and visual clarity. The vector roads were transformed from the raster traffic and classified into levels 1–8 based on their traffic volumes, from high to low: Level 1 (40001 and above), Level 2 (20001–40000), Level 3 (12001–20000), Level 4 (8001–12000), Level 5 (3001–8000), Level 6 (501–3000), Level 7 (101–500), and Level 8 (1–100) (Fig. 4).

Fig. 4

Road network maps during the Bronze and early Iron Ages in the Haidai region of China. a Shang Dynasty; b Western Zhou Dynasty; c Spring & Autumn Period; d Warring States Period

The Shang Dynasty’s road structure was simple, with only short-distance roads. The Shang Dynasty’s main roads extended eastward from the Yin Ruins in Anyang, Central China, bypassing the Taiyi Mountains. They were divided into two eastward, one leading to the Qianzhangda site and the other to the Daxinzhuang and Subutun sites [35]. The roads became more complicated from Western Zhou, Spring & Autumn to the Warring States Period. They developed significantly faster than those in the Shang Dynasty (Fig. 4), with an increase in road length and average traffic volume in some areas such as southwest and southeast Shandong Province, as well as the northern side of the Taiyi Mountains. Finally, dense roads covered almost the entire Haidai region during the Warring States Period. The road network emerged, as evidenced by an increase in road length, complexity of the road structure, and the formation of long-distance communication roads centered on Linzi (the capital of Qi State), Qufu (the capital of Lu State), Ju (the capital of Ju State), and so on. Furthermore, north–south communications were significantly improved through the central Taiyi Mountains, which aligned with findings from archaeological investigations [71].

Table 1 The length and proportion of roads with different levels in the Haidai region of China during the Bronze and Early Iron Ages

The road lengths were calculated by multiplying the number of raster cells representing roads at different levels by the raster’s spatial resolution (Table 1, Fig. 5). The total lengths were 1023,231.87 km, with Shang Dynasty accounting for 232,939.35 km, Western Zhou for 273,453.84 km, Spring & Autumn for 23,129.91 km, and Warring States for 285,547.77 km. However, level 8 roads were excluded due to low traffic, implying little flow. The total road lengths (sum of first-7th level roads) were 19111.86 km for the Shang Dynasty, 35269.02 km for Western Zhou, 51,555.42 km for Spring & Autumn, and 77456.88 km for the Warring States. Furthermore, the proportion is calculated by dividing the length of roads at each level (1st-7th) by the total road length. From the Shang Dynasty to the Warring States Period, the length of each road hierarchy increased almost entirely, with an unusual reduction in first-grade roads, most likely due to the disappearance of coastal routes on the Jiaodong Peninsula. Overall, the changes in road length matched the proportion of the first six levels. However, on the seventh level road from Shang Dynasty to Warring States, the length increased while the proportion decreased, with percentages of 809, 651, 602, and 556‰, respectively.

Fig. 5

Quantitative statistics of roads of various levels in the Haidai region of China during the Bronze and Early Iron Ages. a Road lengths; b Proportion of roads

Fig. 6

Comparison of the simulated road network and road control points belonging to the Warring States Period in the Haidai region of China

Fig. 7

The sum of buffered traffic values extracted from random points (mean value of 41,840) compared to actual point values (200,211) in the flow network (Z-score > 7)

Discussion

Robustness checks on simulated ancient road networks

First, some related relic sites were used to evaluate the accuracy of the SIM-generated road networks, which included the arsenal, ancient battlefield, and Great Wall pass, all of which were likely located at critical road network nodes (Fig. 6). A total of 64 relic sites, hereinafter referred to as control points, were collected from historical records and archaeological excavation reports, including 23 arsenals, 16 ancient battlefields, and 25 Great Wall passes belonging to the Qi and Lu states [72, 73]. Furthermore, 56 control points were preserved because the remaining 8 points overlapped with the site points used to create these ancient road networks. Simple statistics showed that 39 (69.64%) road control points fell on the simulated road network, which is higher than the 57.36% accuracy achieved by Frachetti et al. [12] under similar assessments. There are several possible explanations for this nonconformity. On the one hand, the least cost principle may cause intra-regional traffic to congregate on the least expensive road rather than dispersing across multiple roads [26]. On the other hand, some road control points, such as those along the southeastern coast, are strategically located rather than in densely populated areas.

Further quantitative studies and validation were conducted to assess the model’s statistical accuracy. The statistical results of this test revealed that the mean value of the sum of the simulated traffic was 41,840, the minimum and maximum traffic values were 6585 and 134,621, respectively, and the standard deviation (s.d.) was 21,822 (Fig. 7). In contrast, the traffic value was 200,211 calculated using control points, significantly higher than the one reconstructed using random points. According to the calculated Z score on the random point queue (Z Score = 7.26), the traffic values simulated by the control points far exceeded the mean value calculated by the random points (more than seven standard deviations). The one-sample t-test (p < 0.01) confirmed a significant difference in simulated flow values between road control points and random points. In a nutshell, the high traffic roads simulated by SIM exhibit a non-random good spatial correspondence to historically recorded sites.

Evolution of roads during the bronze and early iron ages in the Haidai region of China

In the Shang Dynasty, the total length of roads (average road flow) was 19,112 km (6.6), and then increased to 35,269 km (31.7) in the Western Zhou Dynasty, 51,555 km (42.8) in the Spring & Autumn Period, and 77,456 km (75.5) in the Warring States Period. Road networks in the Haidai region evolved from relatively simple routes during the Shang Dynasty to a more complex network during the Zhou Dynasty (Figs. 4 and 8).

Fig. 8

The line densities of road networks in the Haidai region of China during the Bronze and Early Iron Ages. a Shang Dynasty; b Western Zhou Dynasty; c Spring & Autumn Period; d Warring States Period

During the Shang Dynasty, roads with rank 7 made up 80.9% of the total length of roads with ranks 1–7, and these roads were likely similar to modern inter-village paths. For example, the Daxinzhuang site was a major road network hub. However, the activity radius was only 5 km because the primary resources, such as stone, were only extracted from the neighboring mountainous area [74], and the artisanal products most likely only met the self-sufficiency needs represented by the bone-making industry [75]. The small-scale economy confined to a short distance made establishing a regional road network difficult.

Nonetheless, there were a few faint long-distance road lines in the Shang Dynasty. As the capital of the Shang dynasty, Yin was one of the largest cities during the Bronze Age. In contrast, the Daxinzhuang site, significantly smaller in scale compared to Yin, was still one of the most prominent in the Haidai region [34] despite not being a city at that time [76]. A solid urban gravitational force of Yin resulted in the formation of a centripetal transportation network with six levels of roads. One road extended northeastward along the Yimeng Mountains, passing through significant sites such as Daxinzhuang and Subutun (Fig. 4a), confirming the historical route of the “Eastern expedition to the barbarians” [77]. The scarcity of produced salt in Yin and its surroundings necessitated extensive transportation of sea salt from northern Shandong Province [78], which most likely prompted the development of these long-distance roads from the eastern coast to Yin [77].

According to historians, the road network in China was formed during the Spring & Autumn Period [79]. Our findings showed that the road networks in the Western Zhou Dynasty and Spring & Autumn Period were similar (Fig. 4b, c). The length of roads at levels 1–6 increased to 34.86 and 39.79% of total road lengths in the Western Zhou Dynasty and the Spring & Autumn Period, respectively, indicating significant transportation progress from the Shang Dynasty to the Spring & Autumn Period. Road networks eventually emerged, focusing on state capitals in the Haidai region such as Linzi, Qufu, and Ju. At the time, this area had more than one-sixth of China’s large cities, with 69 in Lu State and 46 in Qi State [80]. The growth of cities facilitated communication and the formation of radial road networks centered on major cities [81]. Furthermore, the road network broke through topographic barriers, forming roads that crossed mountainous terrain in central Shandong Province [71].

The road network expanded during the Warring States Period, reaching 77,456 km, with the average flow nearly doubling that of the Spring & Autumn Period. The total length of class 1–3 roads increased significantly during the Warring States Period compared to the Shang Dynasty, Western Zhou Dynasty, and Spring & Autumn Period. This rapid development of the road network was most likely attributed to significant social transformations during the Warring States Period. Iron was widely used in agriculture and the military [43, 82, 83], and the popularity of iron farming tools helped significantly increase agricultural productivity [43, 84], promoting food and population growth. The population increased from over 10 million in the Spring & Autumn Period to over 30 million during the Warring States Period [32, 66, 85]. As the population grew, so did the number of cities, which surpassed thousands [37], and this phase saw the second peak of city construction [36]. The development of road networks may had been aided by the growing population and the growth of cities in the Haidai region during the Warring States Period, which benefited from the development of iron smelting, sea salt industry, and commerce. Linzi, Qi State’s capital, was most likely the birthplace of iron smelting [86] and one of several major iron smelting centers in China during the Warring States Period [39]. Furthermore, the 3200 km-long coastlines most likely aided the development of the sea salt industry [78], as 250 sites were discovered in seven major site groups in the lower Xiaoqing River belonging to Qi State [40]. Merchants, represented by “merchant sage” Fan Li, primarily traded salt and iron in this region. [43, 83]. The development of iron smelting, salt production, and commerce all contributed to expanding the road network in this area.

Furthermore, the overall length of road networks increased and communication improved from the Western Zhou Dynasty to the Warring States Period, while the mean exchange distances decreased (Table 2), possibly due to city growth, social unrest, and state annexation. The second peak of city building in the Eastern Zhou Dynasty resulted in the expansion of medium and large settlements in the Haidai region and the emergence of more regional higher-order settlement centers [36]. The proliferation of subordinate higher-order settlement centers most likely resulted in increased intra-regional exchange rather than long-distance exchange. Meanwhile, the East Zhou period was regarded as the most turbulent and warlike in Chinese history [87]. In the early Spring & Autumn Period, there were still over 1200 vassal states, with only seven major states remaining in late Warring States [85]. The great states facilitated domestic exchanges while limiting interstate communication [87, 88]. The emergence of multiple sub-advanced settlement centers and national border limitations are likely to have contributed to the reduction in urban exchange distances.

Table 2 The average interaction radius (km) of representative settlements in the Haidai region of China during the Bronze and Early Iron Ages

Conclusion

The road networks in the Haidai region between the Bronze Age and the Early Iron Age were reconstructed and analyzed using SIM. Regarding location and statistical patterns, the simulated road networks closely matched the road relic points.

Throughout the Bronze Age and Early Iron Age, communication networks in the Haidai area were constantly developing and changing. The road’s total length increased from 19,112 km in the Shang Dynasty to 35,269 km in the Western Zhou Dynasty, 51,555 km in the Spring & Autumn Period, and 77,456 km in the Warring States Period, with the average road traffic of 6.6, 31.7, 42.8, and 75.5, respectively. Roads were still simple during the Shang Dynasty, but the attraction of Yin in the Central Plains prompted the centripetal routes that originated from Yin. Population growth and urbanization aided in the formation of road networks during the Western Zhou Dynasty and the Spring & Autumn Period, while the growth and prosperity of iron smelting, salt-making industry, and commerce fueled the rapid development of road networks during the Warring States Period. However, the gravity of urban exchange distances decreased because of frequent wars and the spread of sub-high-order settlements.

Our study presents the quantitative reconstruction of the ancient transportation road network from a spatial archaeological perspective, offering a novel approach to studying ancient cultural exchange. Although the large temporal and spatial scale employed in this study may have an impact on the accuracy of road network reconstruction, this study demonstrates the potential for studying interactive communications and road evolution in ancient times.