Gravity law in the Chinese highway freight transportation networks

The gravity law has been documented in many socioeconomic networks, which states that the flow between two nodes positively correlates with the strengths of the nodes and negatively correlates with the distance between the two nodes. However, such research on highway freight transportation networks (HFTNs) is rare. We construct the directed and undirected highway freight transportation networks between 338 Chinese cities using about 15.06 million truck transportation records in five months and test the traditional and modified gravity laws using GDP, population, and per capita GDP as the node strength. It is found that the gravity law holds over about two orders of magnitude for the whole sample, as well as the daily samples, except for the days around the Spring Festival during which the daily sample sizes are significantly small. Accordingly, the daily exponents of the gravity law are stable except during the Spring Festival period. The results also show that the gravity law has higher explanatory power for the undirected HFTNs than for the directed HFTNs. However, the traditional and modified gravity laws have comparable explanatory power.

As far as highway network is concerned, Jung, Wang and Stanley investigated the Korean highway network between the 30 largest cities [24]. They found that the traffic between two cities is positively proportional to the populations and negatively proportional to the distance, fitting exactly the traditional gravity law over two orders of magnitude. Kwon and Jung investigated the express bus flow in Korea consisting of 74 cities and 170 bus routes with 6692 operating buses per day and confirmed the presence of the gravity law [25]. Hong and Jung studied the urban bus networks of Korean cities and confirmed the presence of the gravity law [26].
In this work, we investigate a huge database about the freight highway transportation by trucks between 338 cities in mainland China. Although most studies dealt with undirected transportation networks [27,28], we also consider directed transportation networks due to the availability of data. In our analysis, we consider the traditional gravity law, as well as the modified form of the gravity law in which the power-law exponents of the variables are not necessarily fixed.

Data description
The database we analyze was provided by a leading truck logistics company in China. The trucks transport freights between 338 cities in mainland China. The data covers the period from 1 January 2019 to 31 May 2019 and there are totally 15.06 million truck freight transportation records. Each record contains the origin and destination cities and the starting date of the transportation.
We construct the daily flow matrix F(t) = [F ij (t)] 338×338 for each day t, where F ij (t) stands for the number of trucks with freights departing city i for city j on day t. Running trucks that do not load freights are not counted in. The flow matrix for the whole sample is thus We note that F(t) and F are directed networks. When we consider the total truck flow W ij , we obtain an undirected network: where triu denotes the upper triangle matrix operator and T is the transpose operator. By definition, we have In other words, intracity transportation is not included in the data.

The traditional gravity model
In the perspective of directed complex networks, the traditional gravity model states that the edge flow F ij between two nodes i and j takes the form of for i = j, where C is a constant, d ij stands for the distance between i and j, and M i and M j stand respectively for the economic dimensions of the two nodes i and j that are being measured. For undirected transportation networks, we have for j > i. In our analysis, we consider M to be GDP (G), population (P) and per capita GDP (G/P), respectively. Figure 1(a) shows the scatter plot of F ij with respect to G i G j /d 2 ij for the directed network. When M in Eq. (4) stands for GDP, the regression equation of the gravity law reads

Whole transportation networks
where is the error term. We obtain that α 0 = 1.383 ± 0.004 and α 1 = 0.612 ± 0.004, where the adjusted R 2 statistic is 0.476, the F statistic and p-value for the full model are respectively 85,246 and 0.000, and an estimate of the error variance is 0.352. We bin the data with respect to F and illustrate the results in Fig. 1(b). We observe that there is a power-law dependence when F is greater than about 100. A regression shows that α 0 = 1.489±0.025 and α 1 = 0.965 ± 0.017 with their p-values being 0.0000 and 0.0000. The adjusted R 2 statistic is 0.995 and the F statistic and p-value for the full model are respectively 3296 and 0.0000. Figure 1(c) shows the scatter plot of W ij with respect to G i G j /d 2 ij for the undirected network. When M in Eq. (5) stands for GDP, the regression equation of the gravity law reads We obtain that α 0 = 1.705 ± 0.005 and α 1 = 0.662 ± 0.005, where the adjusted R 2 statistic is 0.566, the F statistic and p-value for the full model are respectively 66,530 and 0.000, and an estimate of the error variance is 0.309. We bin the data with respect to W and illustrate the results in Fig. 1(d). We observe that there is a power-law dependence when W is greater than about 200. A regression shows that α 0 = 1.917 ± 0.022 and α 1 = 0.828 ± 0.013 with the p-values being 0.0000 and 0.0000. The adjusted R 2 statistic is 0.996 and the F statistic and p-value for the full model are respectively 3912 and 0.0000. We find that, for the scatter plots, the adjusted R 2 statistic for the undirected HFTN (0.566) is greater than that for the directed HFTN (0.476). This result is visible in Fig. 1, which shows that the scatter plot is thinner for W ij when compared with the one for F ij .

Daily transportation networks
We now test the gravity law with daily directed and undirected freight highway transportation networks. We find that most of the daily networks exhibit the gravity law. As an example, the results for the directed network on 15 January 2019 are illustrated in Fig. 2(a). Regression of Eq. (6) for the scatter data points gives that α 0 = 0.196 ± 0.006 and α 1 = 0.194 ± 0.005, where the adjusted R 2 statistic is 0.191, the F statistic and pvalue for the full model are respectively 6779 and 0.000, and an estimate of the error variance is 0.129. Regression for the binning data points shows that α 0 = -0.257 ± 0.130 and α 1 = 0.824 ± 0.080, where the adjusted R 2 statistic is 0.947, the F statistic and p-value for the full model are respectively 449 and 0.000, and an estimate of the error variance is 0.007. The results for the undirected network on 15 January 2019 are illustrated in Fig. 2(b). Regression of Eq. (7) for the scatter data points gives that α 0 = 0.261 ± 0.007 and α 1 = 0.268 ± 0.006, where the adjusted R 2 statistic is 0.281, the F statistic and pvalue for the full model are respectively 8246 and 0.000, and an estimate of the error variance is 0.142. Regression for the binning data points shows that α 0 = -0.047 ± 0.097 and α 1 = 0.801 ± 0.068, where the adjusted R 2 statistic is 0.959, the F statistic and p-value for the full model are respectively 583 and 0.000, and an estimate of the error variance is 0.006.
However, we find that the daily networks around the Chinese New Year (5 February 2019) do not exhibit the gravity law. The transportation flow decreased significantly during the Spring Festival because most of the truck drivers returned home to gather with their families and most companies were also closed. As an example, the results for the directed network on 4 February 2019 are illustrated in Fig. 2(c). Regression of Eq. (6) for the scatter data points gives that α 0 = 0.167 ± 0.011 and α 1 = 0.014 ± 0.008, where the adjusted R 2 statistic is 0.003, the F statistic and p-value for the full model are respectively 13 and 0.000, and an estimate of the error variance is 0.076. Regression for the binning data points shows that α 0 = 0.585 ± 1.437 and α 1 = 0.121 ± 1.247, where the adjusted R 2 statistic is 0.005, the F statistic and p-value for the full model are respectively 0.047 and 0.833, and an estimate of the error variance is 0.118. The results for the undirected network on 15 January 2019 are illustrated in Fig. 2(d). Regression of Eq. (7) for the scatter data points gives that α 0 = 0.176 ± 0.012 and α 1 = 0.034 ± 0.009, where the adjusted R 2 statistic is 0.014, the F statistic and p-value for the full model are respectively 62 and 0.000, and an estimate of the error variance is 0.085. Regression for the binning data points shows that α 0 = -0.642 ± 0.934 and α 1 = 1.142 ± 0.740, where the adjusted R 2 statistic is 0.485, the F statistic and p-value for the full model are respectively 11 and 0.006, and an estimate of the error variance is 0.063. All the four adjusted R 2 statistics are close to zero, implying that the term G i G j /d 2 ij does not have explanatory power for the transportation flow F ij or W ij and the gravity law is absent.
In Fig. 3, we present the evolution of the exponents α 1 of daily directed and undirected freight highway transportation networks. In each case, the exponent fluctuates roughly around a constant. The cone peak or valley corresponds to the dates around the Spring Festival during which the gravity law does not hold.  Figure 4(a) shows the scatter plot of F ij with respect to P i P j /d 2 ij for the directed network. When M in Eq. (4) stands for population, the regression equation of the gravity law reads

Whole transportation networks
where is the error term. We obtain that α 0 = 2.289 ± 0.007 and α 1 = 0.660 ± 0.005, where the adjusted R 2 statistic is 0.441, the F statistic and p-value for the full model are respectively 74,010 and 0.000, and an estimate of the error variance is 0.376. We bin the data with respect to F and illustrate the results in Fig. 4(b). We observe that there is a power-law dependence when F is greater than about 100. A regression shows that α 0 = 2.986±0.012 and α 1 = 1.081 ± 0.019 with their p-values being 0.0000 and 0.0000. The adjusted R 2 statistic is 0.995 and the F statistic and p-value for the full model are respectively 3316 and 0.0000. Figure 4(c) shows the scatter plot of W ij with respect to P i P j /d 2 ij for the undirected network. When M in Eq. (5) stands for GDP, the regression equation of the gravity law reads lg W ij = α 0 + α 1 lg P i P j /d 2 ij + .
We obtain that α 0 = 2.683 ± 0.009 and α 1 = 0.717 ± 0.006, where the adjusted R 2 statistic is 0.525, the F statistic and p-value for the full model are respectively 56,394 and 0.000, and an estimate of the error variance is 0.338. We bin the data with respect to W and illustrate Figure 4 Testing the traditional gravity law with population for the directed and undirected freight highway transportation networks constructed over five months. (a) Scatter plot of F ij with respect to P i P j /d 2 ij for the directed network. (b) Binning plot of F ij with respect to P i P j /d 2 ij for the directed network. (c) Scatter plot of W ij with respect to P i P j /d 2 ij for the undirected network. (d) Binning plot of W ij with respect to P i P j /d 2 ij for the undirected network the results in Fig. 4(d). We observe that there is a power-law dependence when W is greater than about 200. A regression shows that α 0 = 3.202 ± 0.009 and α 1 = 0.919 ± 0.013 with the p-values being 0.0000 and 0.0000. The adjusted R 2 statistic is 0.997 and the F statistic and p-value for the full model are respectively 4647 and 0.0000.
We find that, for the scatter plots, the adjusted R 2 statistic for the undirected HFTN (0.525) is greater than that for the directed HFTN (0.441). This result is visible in Fig. 4, which shows that the scatter plot is thinner for W ij when compared with the one for F ij .

Daily transportation networks
We now test the gravity law with daily directed and undirected freight highway transportation networks. We find that most of the daily networks exhibit the gravity law. As an example, the results for the directed network on 15 January 2019 are illustrated in Fig. 5(a). Regression of Eq. (8) for the scatter data points gives that α 0 = 0.491 ± 0.005 and α 1 = 0.206 ± 0.005, where the adjusted R 2 statistic is 0.176, the F statistic and pvalue for the full model are respectively 6128 and 0.000, and an estimate of the error variance is 0.131. Regression for the binning data points shows that α 0 = 1.043 ± 0.040 and α 1 = 0.939 ± 0.109, where the adjusted R 2 statistic is 0.927, the F statistic and p-value for the full model are respectively 317 and 0.000, and an estimate of the error variance is 0.010. The results for the undirected network on 15 January 2019 are illustrated in Fig. 5(b). Regression of Eq. (9) for the scatter data points gives that α 0 = 0.663 ± 0.007 and α 1 = 0.283 ± 0.007, where the adjusted R 2 statistic is 0.256, the F statistic and p- value for the full model are respectively 7253 and 0.000, and an estimate of the error variance is 0.147. Regression for the binning data points shows that α 0 = 1.203 ± 0.038 and α 1 = 0.896 ± 0.087, where the adjusted R 2 statistic is 0.947, the F statistic and p-value for the full model are respectively 448 and 0.000, and an estimate of the error variance is 0.007.
However, we find that the daily networks around the Chinese New Year (5 February 2019) do not exhibit the gravity law. The transportation flow decreased significantly during the Spring Festival because most of the truck drivers returned home to gather with their families and most companies were also closed. As an example, the results for the directed network on 4 February 2019 are illustrated in Fig. 5(c). Regression of Eq. (8) for the scatter data points gives that α 0 = 0.182 ± 0.009 and α 1 = 0.001 ± 0.008, where the adjusted R 2 statistic is 0.000, the F statistic and p-value for the full model are respectively 0.072 and 0.789, and an estimate of the error variance is 0.076. Regression for the binning data points shows that α 0 = 0.531 ± 0.645 and α 1 = -0.386 ± 1.222, where the adjusted R 2 statistic is 0.047, the F statistic and p-value for the full model are respectively 0.496 and 0.497, and an estimate of the error variance is 0.113. The results for the undirected network on 15 January 2019 are illustrated in Fig. 5(d). Regression of Eq. (9) for the scatter data points gives that α 0 = 0.222 ± 0.010 and α 1 = 0.022 ± 0.010, where the adjusted R 2 statistic is 0.005, the F statistic and p-value for the full model are respectively 20.5 and 0.000, and an estimate of the error variance is 0.086. Regression for the binning data points shows that α 0 = 1.251 ± 0.540 and α 1 = 1.171 ± 1.271, where the adjusted R 2 statistic is 0.251, the F statistic and p-value for the full model are respectively 4.03 and 0.068, and an estimate of the error variance is 0.092. All the four adjusted R 2 statistics are close to zero, implying ij for the undirected network that the term P i P j /d 2 ij does not have explanatory power for the transportation flow F ij or W ij and the gravity law is absent.
In Fig. 6, we present the evolution of the exponents α 1 of daily directed and undirected freight highway transportation networks. In each case, the exponent fluctuates roughly around a constant. The cone peak or valley corresponds to the dates around the Spring Festival during which the gravity law does not hold. Figure 7(a) shows the scatter plot of F ij with respect to (G i /P i )(G j /P j )/d 2 ij for the directed network. When M in Eq. (4) stands for per capita GDP, the regression equation of the gravity law reads

Whole transportation networks
where is the error term. We obtain that α 0 = 4.772 ± 0.030 and α 1 = 0.672 ± 0.006, where the adjusted R 2 statistic is 0.328, the F statistic and p-value for the full model are respectively 45,778 and 0.000, and an estimate of the error variance is 0.452. We bin the data with respect to F and illustrate the results in Fig. 7(b). We observe that there is a power-law dependence when F is greater than about 100. A regression shows that α 0 = 8.961±0.145 and α 1 = 1.496 ± 0.035 with their p-values being 0.0000 and 0.0000. The adjusted R 2 statistic is 0.990 and the F statistic and p-value for the full model are respectively 1830 and 0.0000.

Figure 7
Testing the traditional gravity law with per capita GDP for the directed and undirected freight highway transportation networks constructed over five months. (a) Scatter plot of F ij with respect to We obtain that α 0 = 5.416 ± 0.040 and α 1 = 0.740 ± 0.008, where the adjusted R 2 statistic is 0.387, the F statistic and p-value for the full model are respectively 32,211 and 0.000, and an estimate of the error variance is 0.437. We bin the data with respect to W and illustrate the results in Fig. 7(d). We observe that there is a power-law dependence when W is greater than about 200. A regression shows that α 0 = 8.215 ± 0.143 and α 1 = 1.257 ± 0.035 with the p-values being 0.0000 and 0.0000. The adjusted R 2 statistic is 0.988 and the F statistic and p-value for the full model are respectively 1266 and 0.0000. We find that, for the scatter plots, the adjusted R 2 statistic for the undirected HFTN (0.387) is greater than that for the directed HFTN (0.328). This result is visible in Fig. 7, which shows that the scatter plot is thinner for W ij when compared with the one for F ij .

Daily transportation networks
We now test the gravity law with daily directed and undirected freight highway transportation networks. We find that most of the daily networks exhibit the gravity law. As an example, the results for the directed network on 15 January 2019 are illustrated in   10) for the scatter data points gives that α 0 = 1.185 ± 0.027 and α 1 = 0.185 ± 0.006, where the adjusted R 2 statistic is 0.117, the F statistic and pvalue for the full model are respectively 3788 and 0.000, and an estimate of the error variance is 0.141. Regression for the binning data points shows that α 0 = 5.933 ± 0.624 and α 1 = 1.227 ± 0.156, where the adjusted R 2 statistic is 0.913, the F statistic and p-value for the full model are respectively 262 and 0.000, and an estimate of the error variance is 0.012. The results for the undirected network on 15 January 2019 are illustrated in Fig. 8(b). Regression of Eq. (11) for the scatter data points gives that α 0 = 1.623 ± 0.035 and α 1 = 0.258 ± 0.008, where the adjusted R 2 statistic is 0.175, the F statistic and pvalue for the full model are respectively 4465 and 0.000, and an estimate of the error variance is 0.163. Regression for the binning data points shows that α 0 = 5.947 ± 0.611 and α 1 = 1.188 ± 0.148, where the adjusted R 2 statistic is 0.917, the F statistic and p-value for the full model are respectively 275 and 0.000, and an estimate of the error variance is 0.012.
However, we find that the daily networks around the Chinese New Year (5 February 2019) do not exhibit the gravity law. The transportation flow decreased significantly during the Spring Festival because most of the truck drivers returned home to gather with their families and most companies were also closed. As an example, the results for the directed network on 4 February 2019 are illustrated in Fig. 8(c). Regression of Eq. (10) for the scatter data points gives that α 0 = 0.206 ± 0.042 and α 1 = 0.005 ± 0.009, where the adjusted R 2 statistic is 0.000, the F statistic and p-value for the full model are respectively 1.300 and 0.254, and an estimate of the error variance is 0.076. Regression for the binning data points shows that α 0 = -1.266 ± 11.353 and α 1 = -0.458 ± 2.616, where the adjusted R 2 statistic for the undirected network. (d) Binning data of W ij with respect to (G i /P i )(G j /P j )/d 2 ij for the undirected network is 0.015, the F statistic and p-value for the full model are respectively 0.153 and 0.704, and an estimate of the error variance is 0.117. The results for the undirected network on 15 January 2019 are illustrated in Fig. 8(d). Regression of Eq. (11) for the scatter data points gives that α 0 = 0.330 ± 0.048 and α 1 = 0.027 ± 0.011, where the adjusted R 2 statistic is 0.006, the F statistic and p-value for the full model are respectively 25.2 and 0.000, and an estimate of the error variance is 0.086. Regression for the binning data points shows that α 0 = 4.802 ± 6.791 and α 1 = 0.936 ± 1.580, where the adjusted R 2 statistic is 0.122, the F statistic and p-value for the full model are respectively 1.67 and 0.221, and an estimate of the error variance is 0.108. All the four adjusted R 2 statistics are close to zero, implying that the term (G i /P i )(G j /P j )/d 2 ij does not have explanatory power for the transportation flow F ij or W ij and the gravity law is absent.
In Fig. 9, we present the evolution of the exponents α 1 of daily directed and undirected freight highway transportation networks. In each case, the exponent fluctuates roughly around a constant. The cone peak or valley corresponds to the dates around the Spring Festival during which the gravity law does not hold.

GDP versus population
It is well known that there is an anomalous power-law scaling between the GDP and population at the city level for many countries [29]: where the scaling exponent κ fluctuates for different countries. Particularly, empirical evidence shows that κ = 1.15 ± 0.08 with R 2 = 0.96 for 295 Chinese cities in 2002 [29]. In Fig. 10, we show the GDP and population for the 338 Chinese cities in 2017 on log-log scales. It shows that there is a power-law scaling between G and P. A linear regression of the following equation gives that κ 0 = 0.391 ± 0.083 and κ = 1.116 ± 0.033, where the corresponding p-values are 0.0000 and 0.0000, the adjusted R 2 statistic is 0.771, the F statistic and p-value for the full model are respectively 1132 and 0.0000.

Exponents of GDP and population
From Eq. (9), we have W ij ∼ (P i P j ) α P . Combining Eq. (12) and these two gravity laws for GDP and population, it is easy to obtain that where α G is α 1 in the gravity model for GDP and α P is α 1 in the gravity model for population. We test the relationship (14) with the exponents of the daily transportation networks. The data for the dates during 1 February 2019 to 9 February 2019 are excluded because their corresponding transportation networks do not exhibit the gravity law. In Fig. 11(a), the exponents α G are those in Fig. 3(a) for the directed networks and Fig. 3(c) for the undirected networks, while the exponents α P are those in Fig. 6(a) for the directed networks and Fig. 6(c) for the undirected networks. In Fig. 11(b), the exponents α G are those in Fig. 3(b) for the directed networks and Fig. 3(d) for the undirected networks, while exponents α P are those in Fig. 6(b) for the directed networks and Fig. 6(d) for the undirected networks. It shows that the data points fall around the line α G κ = α P , clsoe to the diagonal.
If we perform linear regression of the data points for each of the cases in Fig. 11, the four slopes are 0.897, 0.892, 0.889 and 0.923, respectively. They are all less than 1 as predicted by Eq. (14). This deviation is partly due to the fact that the derivation of Eq. (14) does not consider the impact of d 2α 1 ij . Figure 11 Testing the relationship α G κ = α P . The data for the dates during 1 February 2019 to 9 February 2019 are excluded. (a) Exponents α G are those in Fig. 3(a) for the directed networks and Fig. 3(c) for the undirected networks, while exponents α P are those in Fig. 6(a) for the directed networks and Fig. 6(c) for the undirected networks. (b) Exponents α G are those in Fig. 3(b) for the directed networks and Fig. 3(d) for the undirected networks, while exponents α P are those in Fig. 6(b) for the directed networks and Fig. 6(d) for the undirected networks Figure 12 Testing the relationship (κ -1)α G/P = α P . The data for the dates during 1 February 2019 to 9 February 2019 are excluded. (a) Exponents α G/P are those in Fig. 9(a) for the directed networks and Fig. 9(c) for the undirected networks, while exponents α P are those in Fig. 6(a) for the directed networks and Fig. 6(c) for the undirected networks. (b) Exponents α G/P are those in Fig. 9(b) for the directed networks and Fig. 9(d) for the undirected networks, while exponents α P are those in Fig. 6(b) for the directed networks and Fig. 6(d) for the undirected networks

Exponents of per capita GDP and population
Combining Eq. (12) and these two gravity laws for GDP and population, it is easy to obtain that where α G is α 1 in the gravity model for GDP and α P is α 1 in the gravity model for population.
We test the relationship (15) with the exponents of the daily transportation networks. The data for the dates during 1 February 2019 to 9 February 2019 are excluded because their corresponding transportation networks do not exhibit the gravity law. In Fig. 12(a), the exponents α G/P are those in Fig. 9(a) for the directed networks and Fig. 9(c) for the undirected networks, while the exponents α P are those in Fig. 6(a) for the directed net-works and Fig. 6(c) for the undirected networks. In Fig. 12, the exponents α G/P are those in Fig. 9(b) for the directed networks and Fig. 9(d) for the undirected networks, while the exponents α P are those in Fig. 6(b) for the directed networks and Fig. 6(d) for the undirected networks.
It shows that the data points in each plot fall around a straight line. However, they are far away from (κ -1)α G/P = α P . This deviation is partly due to the fact that the derivation of Eq. (15) does not consider the impact of d 2α 1 ij . More importantly, it shows that there are other factors influencing the transportation flow. In other words, the underlying "true model" would be more complicated than the gravity law studied in this work.

The modified gravity model
In the perspective of directed complex networks, the modified gravity model states that the edge flow F ij between two nodes i and j takes the form of for i = j, where C is a constant, d ij stands for the distance between i and j, and M i and M j stand respectively for the economic dimensions of the two nodes i and j that are being measured. For undirected transportation networks, we have for j > i. In our analysis, we consider M to be GDP (G), population (P) and per capita GDP (G/P), respectively.

Whole transportation networks
When M in Eq. (16) stands for GDP, the regression equation of the gravity law reads where is the error term. We obtain that α 0 = 1.171 ± 0.063, α 1 = 0.621 ± 0.008, α 2 = 0.637±0.009, and α 3 = 1.191±0.013, where the adjusted R 2 statistic is 0.477, the F statistic and p-value for the full model are respectively 28,447 and 0.000, and an estimate of the error variance is 0.352. Figure 13(a) shows the scatter plot of F ij with respect to Gα 1 i Gα 2 j dα 3 ij for the directed network. We bin the data with respect to F and illustrate F ij against Gα 1 i Gα 2 j dα 3 ij in Fig. 13(b). We observe that there is a power-law dependence when F is greater than about 30: A regression shows that α 0 = 1.153 ± 0.031 and α 1 = 1.578 ± 0.028 with their p-values being 0.0000 and 0.0000. The adjusted R 2 statistic is 0.995, and the F statistic and p-value for the full model are respectively 3258 and 0.0000. When M in Eq. (17) stands for GDP, the regression equation of the gravity law reads Figure 13 Testing the modified gravity law with GDP for the directed and undirected freight highway transportation networks constructed over five months. (a) Scatter plot of F ij with respect to G We obtain that α 0 = 1.273 ± 0.081, α 1 = 0.743 ± 0.011, α 2 = 0.649 ± 0.011, and α 3 = 1.255 ± 0.016, where the adjusted R 2 statistic is 0.568, the F statistic and p-value for the full model are respectively 22,379 and 0.000, and an estimate of the error variance is 0.308. Figure 13(c) shows the scatter plot of W ij with respect to Gα 1 i Gα 2 j dα 3 ij for the undirected network. We bin the data with respect to F and illustrate F ij against Gα 1 i Gα 2 j dα 3 ij in Fig. 13(d). We observe that there is a power-law dependence when F is greater than about 30: We bin the data with respect to W and illustrate the results in Fig. 13(d). We observe that there is a power-law dependence when W is greater than about 200. A regression shows that α 0 = 1.314 ± 0.030 and α 1 = 1.289 ± 0.021 with the p-values being 0.0000 and 0.0000. The adjusted R 2 statistic is 0.996, and the F statistic and p-value for the full model are respectively 3791 and 0.0000. We find that, for the scatter plots, the adjusted R 2 statistic for the undirected HFTN (0.568) is greater than that for the directed HFTN (0.477). This result is visible in Fig. 13, which shows that the scatter plot is thinner for W ij when compared with the one for F ij .

Daily transportation networks
We now test the modified gravity law with daily directed and undirected freight highway transportation networks. We find that most of the daily networks exhibit the modified  Fig. 14(a). Regression of Eq. (18) for the scatter data points gives that α 0 = -0.119 ± 0.062, α 1 = 0.242 ± 0.010, α 2 = 0.202 ± 0.010, and α 3 = 0.347 ± 0.012, where the adjusted R 2 statistic is 0.195, the F statistic and p-value for the full model are respectively 2315 and 0.000, and an estimate of the error variance is 0.128. Regression of Eq. (19) for the binning data points shows that α 0 = -1.604 ± 0.233 and α 1 = 4.240 ± 0.370, where the adjusted R 2 statistic is 0.957, the F statistic and p-value for the full model are respectively 556 and 0.000, and an estimate of the error variance is 0.006. The results for the undirected network on 15 January 2019 are illustrated in Fig. 14(b). Regression of Eq. (20) for the scatter data points gives that α 0 = -0.139 ± 0.078, α 1 = 0.333 ± 0.013, α 2 = 0.274 ± 0.012, and α 3 = 0.482 ± 0.016, where the adjusted R 2 statistic is 0.286, the F statistic and p-value for the full model are respectively 2815 and 0.000, and an estimate of the error variance is 0.141. Regression of Eq. (21) for the binning data points shows that α 0 = -1.241 ± 0.180 and α 1 = 2.977 ± 0.233, where the adjusted R 2 statistic is 0.965, the F statistic and p-value for the full model are respectively 694 and 0.000, and an estimate of the error variance is 0.005.
However, we find that the daily networks around the Chinese New Year (5 February 2019) do not exhibit the gravity law. The transportation flow decreased significantly during the Spring Festival because most of the truck drivers returned home to gather with their families and most companies were also closed. As an example, the results for the directed network on 4 February 2019 are illustrated in Fig. 14(c). Regression of Eq. (18) for the scatter data points gives that α 0 = -0.342 ± 0.106, α 1 = 0.067 ± 0.017, α 2 = 0.049 ± 0.017, and α 3 = -0.041 ± 0.021, where the adjusted R 2 statistic is 0.022, the F statistic and p-value for the full model are respectively 36.6 and 0.000, and an estimate of the error variance is 0.074. Regression of Eq. (19) for the binning data points shows that α 0 = -11.425 ± 10.233 and α 1 = 22.681 ± 19.102, where the adjusted R 2 statistic is 0.412, the F statistic and pvalue for the full model are respectively 6.999 and 0.024, and an estimate of the error variance is 0.070. The results for the undirected network on 15 January 2019 are illustrated in Fig. 14(d). Regression of Eq. (20) for the scatter data points gives that α 0 = -0.369 ± 0.118, α 1 = 0.085 ± 0.019, α 2 = 0.077 ± 0.019, and α 3 = -0.006 ± 0.023, where the adjusted R 2 statistic is 0.032, the F statistic and p-value for the full model are respectively 48.8 and 0.000, and an estimate of the error variance is 0.084. Regression of Eq. (21) for the binning data points shows that α 0 = -9.981 ± 3.990 and α 1 = 17.909 ± 6.637, where the adjusted R 2 statistic is 0.742, the F statistic and p-value for the full model are respectively 34.563 and 0.000, and an estimate of the error variance is 0.032. The two adjusted R 2 statistics for the scatter data are close to zero, implying that the term G α 1 ij does not have explanatory power for the transportation flow F ij or W ij and the modified gravity law is absent.
In Fig. 15, we present the evolution of the exponents α 1 of daily directed and undirected freight highway transportation networks. In each case, the exponent fluctuates roughly around a constant. The cone peak or valley corresponds to the dates around the Spring Festival during which the gravity law does not hold.

Whole transportation networks
When M in Eq. (16) stands for population, the regression equation of the gravity law reads where is the error term. We obtain that α 0 = 1.409 ± 0.066, α 1 = 0.789 ± 0.011, α 2 = 0.721±0.012, and α 3 = 1.192±0.013, where the adjusted R 2 statistic is 0.446, the F statistic and p-value for the full model are respectively 25,128 and 0.000, and an estimate of the error variance is 0.372. Figure 16(a) shows the scatter plot of F ij with respect to Pα 1 i Pα 2 j dα 3 ij for the directed network. We bin the data with respect to F and illustrate F ij against Pα 1 i Pα 2 j dα 3 ij in Fig. 16(b). We observe that there is a power-law dependence when F is greater than about 30: A regression shows that α 0 = 1.539 ± 0.024 and α 1 = 1.639 ± 0.028 with their p-values being 0.0000 and 0.0000. The adjusted R 2 statistic is 0.995, and the F statistic and p-value for the full model are respectively 3453.6 and 0.0000.
We obtain that α 0 = 1.593 ± 0.085, α 1 = 0.868 ± 0.014, α 2 = 0.797 ± 0.015, and α 3 = 1.271 ± 0.017, where the adjusted R 2 statistic is 0.531, the F statistic and p-value for the full model are respectively 19,295 and 0.000, and an estimate of the error variance is 0.334. Figure 16(c) shows the scatter plot of W ij with respect to Pα 1 i Pα 2 j dα 3 ij for the undirected network. We bin the data with respect to F and illustrate F ij against Pα 1 i Pα 2 j dα 3 ij in Fig. 16(d). We observe that there is a power-law dependence when F is greater than about 30: We bin the data with respect to W and illustrate the results in Fig. 16(d). We observe that there is a power-law dependence when W is greater than about 200. A regression shows that α 0 = 1.741 ± 0.024 and α 1 = 1.333 ± 0.022 with the p-values being 0.0000 and 0.0000. The adjusted R 2 statistic is 0.995, and the F statistic and p-value for the full model are respectively 3631 and 0.0000. We find that, for the scatter plots, the adjusted R 2 statistic for the undirected HFTN (0.531) is greater than that for the directed HFTN (0.446). This result is visible in Fig. 16, which shows that the scatter plot is thinner for W ij when compared with the one for F ij .

Daily transportation networks
We now test the modified gravity law with daily directed and undirected freight highway transportation networks. We find that most of the daily networks exhibit the modified gravity law. As an example, the results for the directed network on 15 January 2019 are illustrated in Fig. 17(a). Regression of Eq. (22) for the scatter data points gives that α 0 = -0.244 ± 0.066, α 1 = 0.353 ± 0.014, α 2 = 0.251 ± 0.014, and α 3 = 0.336 ± 0.013, where the adjusted R 2 statistic is 0.192, the F statistic and p-value for the full model are respectively 2278 and 0.000, and an estimate of the error variance is 0.129. Regression of Eq. (23) for the binning data points shows that α 0 = -2.188 ± 0.319 and α 1 = 4.334 ± 0.426, where the adjusted R 2 statistic is 0.946, the F statistic and p-value for the full model are respectively 439 and 0.000, and an estimate of the error variance is 0.007. The results for the undirected network on 15 January 2019 are illustrated in Fig. 17(b). Regression of Eq. (24) for the scatter data points gives that α 0 = -0.189 ± 0.083, α 1 = 0.420 ± 0.017, α 2 = 0.361 ± 0.017, and α 3 = 0.470 ± 0.016, where the adjusted R 2 statistic is 0.271, the F statistic and p-value for the full model are respectively 2616 and 0.000, and an estimate of the error variance is 0.144. Regression of Eq. (25) for the binning data points shows that α 0 = -1.427 ± 0.210 and α 1 = 3.046 ± 0.257, where the adjusted R 2 statistic is 0.960, the F statistic and p-value for the full model are respectively 596 and 0.000, and an estimate of the error variance is 0.006.
However, we find that the daily networks around the Chinese New Year (5 February 2019) do not exhibit the gravity law. The transportation flow decreased significantly during the Spring Festival because most of the truck drivers returned home to gather with their families and most companies were also closed. As an example, the results for the directed  23) for the binning data points shows that α 0 = -9.526 ± 13.747 and α 1 = 24.918 ± 33.416, where the adjusted R 2 statistic is 0.216, the F statistic and pvalue for the full model are respectively 2.76 and 0.128, and an estimate of the error variance is 0.093. The results for the undirected network on 15 January 2019 are illustrated in Fig. 17(d). Regression of Eq. (24) for the scatter data points gives that α 0 = -0.242 ± 0.128, α 1 = 0.097 ± 0.027, α 2 = 0.066 ± 0.027, and α 3 = -0.005 ± 0.024, where the adjusted R 2 statistic is 0.017, the F statistic and p-value for the full model are respectively 24.7 and 0.000, and an estimate of the error variance is 0.085. Regression of Eq. (25) for the binning data points shows that α 0 = -12.333 ± 5.064 and α 1 = 28.166 ± 10.873, where the adjusted R 2 statistic is 0.726, the F statistic and p-value for the full model are respectively 31.854 and 0.000, and an estimate of the error variance is 0.034. The two adjusted R 2 statistics for the scatter data are close to zero, implying that the term P α 1 ij does not have explanatory power for the transportation flow F ij or W ij and the modified gravity law is absent.
In Fig. 18, we present the evolution of the exponents α 1 of daily directed and undirected freight highway transportation networks. In each case, the exponent fluctuates roughly around a constant. The cone peak or valley corresponds to the dates around the Spring Festival during which the gravity law does not hold.
ij for the undirected network

Whole transportation networks
When M in Eq. (16) stands for per capita GDP, the regression equation of the modified gravity law reads where is the error term. We obtain that α 0 = 4.984 ± 0.050, α 1 = 0.544 ± 0.018, α 2 = 0.681 ± 0.019, and α 3 = 1.385 ± 0.014, where the adjusted R 2 statistic is 0.330, the F statistic and p-value for the full model are respectively 15,365 and 0.000, and an estimate of the error variance is 0.451. Figure 19(a) shows the scatter plot of F ij with respect to (G i /P i )α 1 (G j /P j )α 2 /dα 3 ij for the directed network. We bin the data with respect to F and illustrate F ij against (G i /P i )α 1 (G j /P j )α 2 /dα 3 ij in Fig. 19(b). We observe that there is a power-law dependence when F is greater than about 30: A regression shows that α 0 = 9.317 ± 0.158 and α 1 = 2.191 ± 0.053 with their p-values being 0.0000 and 0.0000. The adjusted R 2 statistic is 0.990, and the F statistic and p-value for the full model are respectively 1736 and 0.0000.

Figure 19
Testing the modified gravity law with population for the directed and undirected freight highway transportation networks constructed over five months. (a) Scatter plot of F ij with respect to ij for the directed network. (c) Scatter plot of W ij with respect to (G i /P i ) α 1 (G j /P j ) α 2 /d We obtain that α 0 = 5.621 ± 0.067, α 1 = 0.660 ± 0.026, α 2 = 0.701 ± 0.023, and α 3 = 1.520 ± 0.019, where the adjusted R 2 statistic is 0.388, the F statistic and p-value for the full model are respectively 10,769 and 0.000, and an estimate of the error variance is 0.436. Figure 19(c) shows the scatter plot of W ij with respect to ( for the undirected network. We bin the data with respect to F and illustrate F ij against (G i /P i )α 1 (G j /P j )α 2 /dα 3 ij in Fig. 19(d). We observe that there is a power-law dependence when F is greater than about 30: We bin the data with respect to W and illustrate the results in Fig. 19(d). We observe that there is a power-law dependence when W is greater than about 200. A regression shows that α 0 = 8.641 ± 0.140 and α 1 = 1.726 ± 0.043 with the p-values being 0.0000 and 0.0000. The adjusted R 2 statistic is 0.990, and the F statistic and p-value for the full model are respectively 1642 and 0.0000.
We find that, for the scatter plots, the adjusted R 2 statistic for the undirected HFTN (0.388) is greater than that for the directed HFTN (0.330). This result is visible in Fig. 19, which shows that the scatter plot is thinner for W ij when compared with the one for F ij .

Daily transportation networks
We now test the modified gravity law with daily directed and undirected freight highway transportation networks. We find that most of the daily networks exhibit the modified gravity law. As an example, the results for the directed network on 15 January 2019 are illustrated in Fig. 20(a). Regression of Eq. (26) for the scatter data points gives that α 0 = 1.147 ± 0.043, α 1 = 0.199 ± 0.019, α 2 = 0.197 ± 0.019, and α 3 = 0.364 ± 0.013, where the adjusted R 2 statistic is 0.117, the F statistic and p-value for the full model are respectively 1264 and 0.000, and an estimate of the error variance is 0.141. Regression of Eq. (27) for the binning data points shows that α 0 = 5.689 ± 0.583 and α 1 = 6.644 ± 0.830, where the adjusted R 2 statistic is 0.916, the F statistic and p-value for the full model are respectively 272 and 0.000, and an estimate of the error variance is 0.012. The results for the undirected network on 15 January 2019 are illustrated in Fig. 20(b). Regression of Eq. (28) for the scatter data points gives that α 0 = 1.582 ± 0.056, α 1 = 0.271 ± 0.024, α 2 = 0.272 ± 0.023, and α 3 = 0.509 ± 0.017, where the adjusted R 2 statistic is 0.175, the F statistic and p-value for the full model are respectively 1489 and 0.000, and an estimate of the error variance is 0.163. Regression of Eq. (29) for the binning data points shows that α 0 = 5.773 ± 0.584 and α 1 = 4.622 ± 0.569, where the adjusted R 2 statistic is 0.918, the F statistic and p-value for the full model are respectively 280 and 0.000, and an estimate of the error variance is 0.011. However, we find that the daily networks around the Chinese New Year (5 February 2019) do not exhibit the gravity law. The transportation flow decreased significantly during the Spring Festival because most of the truck drivers returned home to gather with their families and most companies were also closed. As an example, the results for the directed network on 4 February 2019 are illustrated in Fig. 20(c). Regression of Eq. (26) for the scatter data points gives that α 0 = -0.103 ± 0.071, α 1 = 0.097 ± 0.032, α 2 = 0.124 ± 0.032, and α 3 = -0.039 ± 0.021, where the adjusted R 2 statistic is 0.023, the F statistic and p-value for the full model are respectively 37.4 and 0.000, and an estimate of the error variance is 0.074. Regression of Eq. (27) for the binning data points shows that α 0 = -8.942 ± 4.624 and α 1 = 32.579 ± 15.580, where the adjusted R 2 statistic is 0.685, the F statistic and pvalue for the full model are respectively 21.7 and 0.001, and an estimate of the error variance is 0.037. The results for the undirected network on 15 January 2019 are illustrated in Fig. 20(d). Regression of Eq. (28) for the scatter data points gives that α 0 = 0.007 ± 0.079, α 1 = 0.122 ± 0.037, α 2 = 0.156 ± 0.034, and α 3 = 0.004 ± 0.023, where the adjusted R 2 statistic is 0.028, the F statistic and p-value for the full model are respectively 42.6 and 0.000, and an estimate of the error variance is 0.084. Regression of Eq. (29) for the binning data points shows that α 0 = -2.627 ± 2.288 and α 1 = 15.534 ± 10.407, where the adjusted R 2 statistic is 0.468, the F statistic and p-value for the full model are respectively 10.6 and 0.007, and an estimate of the error variance is 0.065. The two adjusted R 2 statistics for the scatter data are close to zero, implying that the term (G i /P i ) α 1 (G j /P j ) α 2 /d α 3 ij does not have explanatory power for the transportation flow F ij or W ij and the modified gravity law is absent.
In Fig. 21, we present the evolution of the exponents α 1 of daily directed and undirected freight highway transportation networks. In each case, the exponent fluctuates roughly around a constant. The cone peak or valley corresponds to the dates around the Spring Festival during which the gravity law does not hold.

Summary and discussions
In this work, we constructed the directed and undirected highway freight transportation networks between 338 Chinese cities using about 15.06 million truck transportation records in five months (2019/01/05-2019/05/31) and tested the traditional and modified gravity laws using GDP, population and per capita GDP as the node strength. It is found that the gravity law holds over about two orders of magnitude for the whole sample, as well as the daily samples, except for the days around the Spring Festival during which the daily sample sizes are significantly small. Accordingly, the daily exponents of the gravity law are stable except during the Spring Festival period. Other than the Spring Festival holiday, there are a number of holidays in China. However, we observe that there show no abnormal phenomena in the time series of α 1 other than the Spring Festival. These seemingly counterintuitive results are explained by the fact that the Spring Festival is special other than other holidays. It is the tradition that almost all Chinese stayed in other places will return to their hometowns to have dinner with all family members especially parents on the Eve of the Spring Festival. Accordingly, most companies and factories are closed during the Spring Festival for one or two weeks, resulting in a sharp decrease in the demand of freight transportation.
We confirmed that there is a power-law scaling between GDP and population at the city level. We carried out a scaling analysis and obtained the relationship between the expo- ij for the directed network. (b) Binning data of F ij with respect to (G i /P i ) α 1 (G j /P j ) α 2 /d α 3 ij for the directed network. (c) Scatter data of W ij with respect to (G i /P i ) α 1 (G j /P j ) α 2 /d α 3 ij for the undirected network. (d) Binning data of W ij with respect to (G i /P i ) α 1 (G j /P j ) α 2 /d α 3 ij for the undirected network nents of different gravity laws. Using daily transportation networks, the scaling relationship between the GDP-based gravity law and the population-based gravity law has been well verified with mild deviations. The feasible regions for the gravity laws differ from each other (see the figures). However, the scaling relationship could not be confirmed when the per capita GDP-based gravity law is considered. To some extent, this observation is consistent with the fact that the per capita GDP-based gravity law gives significantly smaller adjusted R 2 values.
The adjusted R 2 values show that the gravity law has higher explanatory power for the undirected highway freight transportation networks than for the directed highway freight transportation networks. However, the traditional and modified gravity laws have comparable explanatory power according to the fact that their corresponding adjusted R 2 values are very close to each other. In addition, although the GDP-based and the population based models already have quite high explanatory power for the transportation flows (more than 40% for directed networks and more than 50% for undirected networks), additional factors that account for the unexplained components call for further research. One possible direction is to consider the comparative advantage of different cities embedded in the industrial structure of China as well as different demands from different cities.
We note that the fitted parameters or exponents are not general. The literature has shown that the exponents could be different due to different systems (air lines, trains, buses, et al.), different time periods, and different places. In addition, there are other