Spatial Econometric Analysis of the Relationship between Urban Land and Regional Economic Development in the Beijing–Tianjin–Hebei Coordinated Development Region

: Against the background of coordinated development of the Beijing–Tianjin–Hebei region, it is of great signiﬁcance to quantitatively reveal the contribution rate of the inﬂuencing factors of urban land for optimizing the layout of urban land across regions and innovating the inter-regional urban land supply linkage. However, the interaction e ﬀ ects and spatial e ﬀ ects decomposition have not been well investigated in the existing research studies on this topic. In this study, based on the cross-sectional data in 2015 and using the spatial lag model, spatial error model and spatial Durbin model, we analyzed the relationship between urban land and regional economic development at the county level in the Beijing–Tianjin–Hebei region. The results show that: (1) there are endogenous interaction e ﬀ ects of urban land, and the growth of urban land in a county will drive the corresponding growth of urban land in neighboring counties; (2) the local population, average wages, highway mileage density, and actual utilization of foreign capital have positive e ﬀ ects on the scale of urban land in local and neighboring counties; local GDP in the secondary / tertiary sector and the urbanization rate have positive e ﬀ ects on local urban land scale, but negative e ﬀ ects on the urban land scale of neighboring counties; (3) the contribution degree of the direct e ﬀ ect is ranked as follows: GDP in the secondary / tertiary sector > total population > urbanization rate. The order of factors with a signiﬁcant spatial spillover e ﬀ ect on the scale of urban land in neighboring counties is as follows: average wages > total population > highway mileage density. The GDP in secondary / tertiary sector, population, and urbanization rate are the main inﬂuencing factors for the scale of urban land at the county level in the Beijing–Tianjin–Hebei region. It is an important ﬁnding that average wages are the most prominent among the spatial spillovers. We should attach importance to the spillover e ﬀ ect of geographic space and construct an urban spatial pattern coordinated with economic development.


Introduction
The city is a large-scale residential area formed by the agglomeration of socioeconomic activities, dominated by non-agricultural industries and the agglomeration of non-agricultural population. In November 2014, the Chinese State Council defined a new standard for classifying cities. Cities are divided into super-mega, mega, large, medium, and small cities according to the size of urban permanent residents with a population greater than 10 × 10 6 , between 5 × 10 6 and 10 × 10 6 , between 1 × 10 6 and 5 × 10 6 , between 0.5 × 10 6 and 1 × 10 6 , and below 0.5 × 10 6 , respectively. Urbanization is an important symbol of the level of economic development and reflects the process of social structure reform in a country or region. China is not only the engine of global economic growth but also the engine of global urban population growth. Since the reform and opening up, China's economy has maintained rapid growth for nearly 40 years [1,2], and China's urbanization process has been advancing at an unprecedented speed. The urbanization rate has increased from 17.92% in 1978 to 56.10% in 2015 [3], which means that 771 million people live in urban areas. At the same time, rapid urbanization has accelerated the urbanization of land [4]. Cities attract a large number of immigrants from rural areas due to the higher income and increased job opportunities [5,6]. These migrations, and the growth of the urban population itself, led to a dramatic increase in the urban population. The increase in the urban population and its desire for better housing will inevitably lead to the expansion of urban land. That is, the spread of urban built-up areas in space [7][8][9]. In addition, the rapidly growing urban population increases the urban population density [10]. In the construction of residential areas, it is necessary to pay attention to the surrounding supporting facilities, such as the construction of parks. The land occupied by these supporting facilities is also an important reason for the expansion of urban land. Urban expansion provides land for urban development and effectively supports socioeconomic development [11]. However, the essence of urban land expansion is the conversion of rural and natural land to urban land. This conversion is the transformation of natural land into semi-natural, semi-artificial or artificial land [12][13][14]. Significant changes in the constituent elements of the ecosystem and their quantity-to-ratio relationship cause the original natural surface to lose its original ecological functions, and the resulting ecological environment problems will become one of the biggest problems of the 21st century [15][16][17]; for example, the loss of natural vegetation and farmland [18,19], air pollution [20], water pollution [21], water cycle changes [22], local climate change [23], biodiversity reduction [24], wetland destruction [25], impacts on public health [26], greenhouse gas emissions, and heat island effect, etc. [27,28]. How do we supply the land needed for urban development and minimize the adverse effects of land-use change? Harmonizing the contradiction between the supply and demand of land is the current problem to be considered. Therefore, it is of extraordinary significance to obtain the boundary of urban built-up areas and analyze the influencing factors of the urban land scale to protect limited cultivated land resources and promote the sustainable development of the urban ecological environment [29].
Population growth and economic development are the main factors affecting urban land growth [30,31]. The rapid expansion of urban land in China has been driven not only by economic development and market forces, but also by government plans for urbanization [32]. Urbanization means that the proportion of urban population in the total population increases, which inevitably leads to the increase of land occupied by the urban population. China's rapid urbanization and urban development constantly consume a large number of land resources. The essence of urban land scale growth is to rely on local investment or foreign capital to convert non-agricultural land into construction land. Roads are the lifeblood that connects passenger and freight traffic to different cities. In China, if cities want to develop, transportation must go first. Terrain is the background of urban development, which restricts the growth of urban land. What is the mechanism of these socioeconomic factors affecting the urban land scale? This is the core problem to be solved in this paper.
In China, the administrative region at the county level is the most effective executor of national policies and plays a pivotal role in connecting the above and the following in the administrative region system at the central, provincial, prefecture-level city, county, and township (town) levels. In 2015, the Chinese government issued the Outline of the Beijing-Tianjin-Hebei Coordinated Development Plan, which pointed out that promoting the coordinated development of the Beijing-Tianjin-Hebei is a major national strategy. Taking the Beijing-Tianjin-Hebei as a whole, the adjustment and optimization of urban layout is an important part of coordinated development. Against this background, it is of certain theoretical and practical significance to investigate the influencing factors of urban land. In order to capture the influencing factors of the urban land scale in the Beijing-Tianjin-Hebei Coordinated Development Region (BTHCDR) at a specific time node, this paper adopts the cross-sectional spatial econometrics model to study the influencing factors of the urban land scale in 175 spatial units in Beijing-Tianjin-Hebei in 2015. Variable statistical data are adopted at the county level, including municipal districts, county-level cities and counties, and the central urban areas of prefecture-level cities are not subdivided into municipal districts.
Based on the cross-sectional data of the BTHCDR at the county level in 2015, this paper explores the influence mechanism of socioeconomic factors on the scale of urban land. This paper has three objectives: (1) to investigate whether there is an endogenous interaction effect between urban land scale, an exogenous interaction effect between socioeconomic factors, or an interaction effect between error terms; (2) to quantitatively analyze the influence of each explanatory variable on urban land scale through spatial effects decomposition; (3) to rank the contribution degree of the influencing factors according to the significance and value of the elasticity coefficient, provide valuable insights for the practice of coordinating the urban land layout between regions.
This paper is organized into six sections: Section 1 presents the introduction and research objectives; Section 2, literature review and the actual stage of research on this topic; Section 3 introduces the study area, data source, and research methods; Section 4, analysis of empirical results; Section 5, discussion; Section 6, conclusions, followed by references.

Literature Review
Many literature studies are devoted to exploring the influencing factors of urban land expansion, including natural, socioeconomic, infrastructure, and policy factors. Natural factors are mainly topographic conditions, such as elevation and slope [21,24,28,[33][34][35], that restrict urban land development. Socioeconomic factors include population [21,[34][35][36][37][38], GDP [34,[36][37][38][39], urbanization rate [40,41], fixed-asset investment [34,38,42], foreign direct investment [36,43], and urban wages or disposable income [33,34,36,39,42]. Infrastructure factors are mainly traffic conditions [38,42]. Policy factors are national, provincial, or municipal land-use policies or regulations [28,34,38]. The research methods of the driving force of urban land expansion mainly include multiple linear regression [44][45][46], logistic regression [47][48][49], analytic hierarchy process (AHP) [50], principal component analysis (PCA) [51], ridge regression [38], artificial neural network [52], and cellular automata [24]. The dependent variable of logistic regression is two classified variables, and the independent variable can be a continuous variable or a classified variable, which can reflect the spatial characteristics of the variable. This means where the urban expansion will take place. So it can be used for the driving force analysis and prediction of land-use change. However, the logistic regression model ignores the interaction and random disturbance in the neighborhood. Multiple linear regression uses time series data to explore the relative importance of different factors. Only continuous variables can be used, and it is not applicable when the dependent variable is a classified variable. The analytic hierarchy process (AHP) determines the relative importance of factors at all levels through pairwise comparison, which is a subjective method. The accuracy of the result is largely dependent on expert knowledge, and it is difficult to analyze the spatial information of urban land expansion. The artificial neural network (Ann) has strong nonlinear processing ability. It uses the trained image template to identify the image, and the image extraction results depend on the expert knowledge and experience stored in the image template. Cellular automata (CA) model simulates the process or future scenario of urban expansion by using cells and transformation rules that are discrete in space and time [53].
Many studies have investigated the influencing factors of urban land or construction land scale. However, there are some deficiencies in previous studies. Firstly, the spatial autocorrelation of Sustainability 2020, 12, 8451 4 of 21 spatial units is seldom considered in the study of driving factors of urban expansion. Spatial effects are very important in geographical processes, and ignoring them may lead to biased results [21]. The expansion of urban land in a certain administrative region is not carried out alone. The expansion of urban land is affected by neighboring administrative regions. For example, Zhang et al. [25] compared the spatial and temporal pattern of urban growth in Beijing, Tianjin, and Tangshan from the 1970s to 2013, and emphasized the influence of neighboring cities on the direction of urban expansion. In view of this, scholars suggested that the teleconnections framework of urban land promote urbanization [54]. That is to say, in the process of studying urban land expansion, the spatial autocorrelation of adjacent administrative areas and the layout of regional urban land should be fully considered. Secondly, previous quantitative studies on the driving force of urban land expansion have mostly focused on individual cities or prefecture-level cities as a unit by which to study the national level, provincial level, or the Yangtze River Delta and other economic zones. The explanatory ability of the driving factors determined by a single urban scale to the dependent variables is often difficult to extend to the regional research composed of the county level. The scale of the analysis unit influences the correlation between the land-use system and its potential explanatory variables [55]. The driving effect of explanatory variables on urban land may vary with a change in space and time scale [56]. For example, Li, Li, and Wu [15], based on data from four periods from 1990 to 2008, found that the relative importance of key drivers of urban growth in the central part of the Yangtze River Delta is different at the regional-, prefecture-, and county levels. Thirdly, previous studies mostly analyzed multi-year data and seldom used cross-sectional data, which can capture information of specific time nodes. Gao et al. [57] used panel data from 2000 to 2010 to find that six factors, such as foreign direct investment and population, drove the urban land expansion in the Yangtze River Delta region. Using panel data from 1995 to 2010, Jiang and Zhang [58] found that urban wages were the basic factor to promote the conversion of agricultural land into urban land. Using 2008 cross-sectional data, Ye, Zhang, and Pu [40] concluded that the GDP growth of a region would not only drive the expansion of its own construction land, but also the growth of construction land in neighboring regions. Cross-sectional data have some limitations compared to panel data. For example, the problem of missing variables cannot be solved and the dynamic behavior information of individuals cannot be provided. In addition, because the cross-sectional data have only one dimension, the accuracy of the estimation is not as good as the panel data. However, cross-sectional data can capture specific time node information. In the field of land use, the cross-sectional econometrics model can be used to accurately explore the factors affecting the urban land scale in a certain year. Fourthly, most of the previous studies only analyzed the influence of socioeconomic factors on urban land use or construction land, rarely involving spatial effect decomposition. As for the spatial spillover effect of the influencing factors of urban land or construction land, Wang and Gao have made in-depth research and put forward valuable opinions. Wang [41] analyzed the impact of GDP, population, and urbanization rates on construction land using data from 30 provincial-level sources in mainland China. Based on the 285 Chinese prefecture-level and above cities, Feng et al. [59] investigated the driving forces of urban expansion of different sizes at the national and regional levels. However, in the literature on the topic of influencing factors of urban land or construction land, existing studies have analyzed indirect effect (spatial spillover effect) based on the spatial lag coefficient of explanatory variables. There are few empirical studies that use direct effect, indirect effect, and total effect to deeply analyze the influence of explanatory variables on urban land scale. This research aims to fill this gap.

The Study Area
The study area includes the Beijing, Tianjin, and Hebei provinces (Figure 1), with a total land area of 21.6 × 10 4 km 2 , accounting for about 2.25% of China's total land area. The region has a typical temperate continental monsoon climate: a high temperature and rainy summer, with a cold and dry winter. This area's terrain is northwest high, southeast low, with the Bashang plateau and Yanshan Mountains in the north, and the Taihang Mountains in the west. The east of the Taihang Mountain range and the south of the Yanshan Mountains are part of the north China plain. From west to east, the piedmont alluvial-diluvial plain, mid-lake plain, and coastal plain are located. the piedmont alluvial-diluvial plain, mid-lake plain, and coastal plain are located.
By the end of 2015, the total population of the study area was 111 million, accounting for 8.1% of the total population of China. The GDP was 6.94 trillion yuan, accounting for 10.1% of China's total GDP. The urbanization rate of Beijing, Tianjin, and Hebei was 86.5%, 82.6%, and 51.3%, respectively. The urbanization rate of Hebei is far lower than that of Beijing and Tianjin, and 4.8 percentage points lower than the national average.
The study area is the political and cultural center, and population and economic concentration of China, which has an important strategic position in our country's development. Since 2015, in order to coordinate the development of the Beijing-Tianjin-Hebei city cluster and solve the environmental problems arising from the urbanization process, the state and local governments have issued a number of policies [60].

The Data Source
Urban land is the built-up area of large, medium, and small cities and counties. Built-up area refers to the area which is actually developed and constructed, which is concentrated and contiguous and has basic municipal public facilities. The data of urban land in this paper are from China's 1:100,000 scale remote sensing monitoring database of land use in 2015. The database was independently constructed by Zengxiang Zhang and his team from the land resources remote sensing research department, Aerospace Information Research Institute, Chinese Academy of Sciences. Based on data from Landsat 8 OLI and environment 1 hj-1a CCD, the database was updated based on China's land-use data in 2010.
The space reference adopts the equal area cutting conical projection of Albers. The latitude of the south standard is 25 °N and the latitude of the north standard is 47 °N, the central longitude is 105 °E, the origin of the coordinates is the intersection of 105 °E and the equator, the latitudinal longitude is 0°, and the ellipsoid is Krasovsky [8]. The minimum standard of the figure above is 6 × 6 pixels, and the spatial resolution of 30 m remote sensing image is equivalent to the actual area of 200 × 200 m 2 . The positioning error of the urban land boundary should be less than 0.6 mm, which is equivalent to two pixels in the image or the actual 60 m [61].
Digital elevation data (DEM) are derived from the cold and arid regions science data center (http://westdc.westgis.ac.cn) and used to calculate the topographic relief. Socioeconomic data were collected from the 2016 Hebei Economic Yearbook, Beijing Regional Statistical Yearbook, Beijing Statistical Yearbook and Tianjin Statistical Yearbook. By the end of 2015, the total population of the study area was 111 million, accounting for 8.1% of the total population of China. The GDP was 6.94 trillion yuan, accounting for 10.1% of China's total GDP. The urbanization rate of Beijing, Tianjin, and Hebei was 86.5%, 82.6%, and 51.3%, respectively. The urbanization rate of Hebei is far lower than that of Beijing and Tianjin, and 4.8 percentage points lower than the national average.
The study area is the political and cultural center, and population and economic concentration of China, which has an important strategic position in our country's development. Since 2015, in order to coordinate the development of the Beijing-Tianjin-Hebei city cluster and solve the environmental problems arising from the urbanization process, the state and local governments have issued a number of policies [60].

The Data Source
Urban land is the built-up area of large, medium, and small cities and counties. Built-up area refers to the area which is actually developed and constructed, which is concentrated and contiguous and has basic municipal public facilities. The data of urban land in this paper are from China's 1:100,000 scale remote sensing monitoring database of land use in 2015. The database was independently constructed by Zengxiang Zhang and his team from the land resources remote sensing research department, Aerospace Information Research Institute, Chinese Academy of Sciences. Based on data from Landsat 8 OLI and environment 1 hj-1a CCD, the database was updated based on China's land-use data in 2010.
The space reference adopts the equal area cutting conical projection of Albers. The latitude of the south standard is 25 • N and the latitude of the north standard is 47 • N, the central longitude is 105 • E, the origin of the coordinates is the intersection of 105 • E and the equator, the latitudinal longitude is 0 • , and the ellipsoid is Krasovsky [8]. The minimum standard of the figure above is 6 × 6 pixels, and the spatial resolution of 30 m remote sensing image is equivalent to the actual area of 200 × 200 m 2 . The positioning error of the urban land boundary should be less than 0.6 mm, which is equivalent to two pixels in the image or the actual 60 m [61].
Digital elevation data (DEM) are derived from the cold and arid regions science data center (http://westdc.westgis.ac.cn) and used to calculate the topographic relief. Socioeconomic data were Sustainability 2020, 12, 8451 6 of 21 collected from the 2016 Hebei Economic Yearbook, Beijing Regional Statistical Yearbook, Beijing Statistical Yearbook and Tianjin Statistical Yearbook.

Research Methods
(1) Spatial Autocorrelation Classical correlation analysis only considers the correlation between variables affected by attribute values, while spatial autocorrelation analysis is the correlation between the variables affected by attribute values and spatial positions. Spatial autocorrelation is whether the value of the variable attribute of a spatial unit is related to the value of the same variable attribute of a neighboring spatial unit. In this paper, Moran's I is used to measure the global spatial autocorrelation of variables in the study area. The formula is as follows [62]: where: I is the Moran Index; X i and X j represent the attribute values of space unit i and j, respectively; X is the average value of the attributes of n spatial units; W ij is a spatial weight matrix. If the spatial weight matrix is normalized by row, the value interval of Moran's I is [−1, 1]. Moran's I greater than 0 indicates a positive correlation, and the closer its value is to 1, the more obvious the clustering characteristics of geographic objects. Moran's I less than 0 indicates a negative correlation, and the closer its value is to −1, the more obvious the discrete characteristics of geographic objects are. Moran's I equal to 0 indicates that the geographic objects are randomly distributed in space [63,64]. p value and Z score are commonly used in the significance test of Moran's index. p value tests whether the difference in attribute value is significant, and Z score tests whether the spatial position pattern is significant. With 90% confidence, p value < 0.10, Z score < −1.65 or > +1. 65. With 95% confidence, p < 0.05, Z score < −1.96 or > +1.96. With 99% confidence, p value < 0.01, Z score < −2.58 or > +2. 58. The calculation formula of Z score is as follows: where: E(I) is the expectation of I, and VAR(I) is the variance of I. In this paper, the binary adjacency matrix is established by using Rook adjacency. If there are common edges in the administrative boundary of two counties, the spatial weight is 1, otherwise it is 0, and the spatial weight matrix is row-standardized. In addition to the adjacency matrix, the spatial weight matrix can also be determined according to the distance between each space element. One way is to determine whether two units are adjacent according to the threshold distance, such as inverse distance. Another way is to set the nearest-K units as neighbors, and each unit has K neighbors. For the sake of robustness, the spatial weight matrix based on inverse Euclidean distance is also used to check the reliability of the calculation results.
(2) Spatial Econometric Model The linear regression model is based on the assumption that the error terms are independent of each other and the mean value of the error terms is zero. This assumption is not true when variables are spatially autocorrelated. In this case, the spatial econometric model will be more applicable. Classical spatial econometric models include spatial lag model (SLM), spatial error model (SEM) and spatial Durbin model (SDM) [65].
When the explained variables of a space unit depend on the explained variables of other space units, it is the spatial lag model (SLM), and its formula is: where Y represents the attribute value of the explained variable. If there are n space units in total, Y is the column vector of n × 1. X represents the attribute value of the explanatory variable. If there is k of them, X is an n × k matrix. WY is the spatial lag vector of the explained variable, and if the coefficient ρ is significant, it means that there is an obvious interaction between the explained variables. β is a vector of k × 1, representing the regression coefficient of the explained variable. W represents the spatial weight matrix of (n × n). The ε represents the random error term of n × 1, which satisfies the conditions E(ε) = 0 and cov(ε) = σ 2 I n .
When the error term has a spatial autocorrelation, it is the spatial error model (SEM), the autoregression of the error term, that is considered in the model, and its formula is: where λ is the regression coefficient of the error term. If λ is significant, it indicates that there is indeed a spatial autocorrelation between the error terms. Wµ is the spatial lag vector of the error term. The spatial Durbin model (SDM) takes the exogenous interaction effect of the independent variable and the endogenous interaction effect of the dependent variable into account, namely, the spatial lag term WX of the explanatory variable and the spatial lag term WY of the explained variable [66], whose formula is: where θ is the exogenous interaction effect coefficient of the explanatory variable, whose value reflects the spatial correlation degree of the explanatory variable. The estimation technique used in the linear regression model is the Ordinary Least Squares (OLS) method. The estimation technique used for point estimation in spatial econometric models is the maximum likelihood method. The spatial effects decomposition of SDM adopt the estimation technique proposed by LeSage and Pace in 2009 [67].
(3) Potential Explanatory Variables By referring to relevant research results and considering the availability of data at the county level, factors influencing the size of urban land were selected. In the past, most studies chose GDP to represent economic influencing factors. Since the primary industry is mostly distributed in rural areas, and the secondary and tertiary industries are mostly concentrated in urban areas, this paper chose the secondary and tertiary GDP to represent the economic influencing factors. Population is the main body of all socioeconomic activities and the ultimate service object. As an important content of urban planning, urbanization rate is closely related to urban land scale. Whether urban construction can proceed smoothly depends on investment in fixed assets. Wages determine the level of the population's consumption of education, commerce, housing, and so on. Highway mileage density is a good index to quantitatively reflect traffic conditions. Actual utilization of foreign capital can truly reflect the level of regional foreign capital utilization. Therefore, this paper uses actual utilization of foreign capital to express foreign investment. Previous studies used foreign direct investment (FDI) to represent foreign investment. Different topographic relief has great influence on the difficulty of urban construction. The author expects to quantitatively reveal the influence of these variables on urban land scale. The explained variable is the urban land scale (hectare) in 2015, and the explanatory variables are selected as shown in Table 1.
In Table 1, correlation coefficient refers to the Pearson correlation test results of each explanatory variable and urban land scale. It can be seen that, except for the significant negative correlation of topographic relief, other variables are significantly positive correlated with urban land scale. The single-variable Moran's I of each explanatory variable is listed at the same time, so as to analyze the trend of agglomerated, discrete, or random distribution of variables in space. The spatial distribution of other explanatory variables showed obvious clustering characteristics, except for the total population and highway mileage density.

Identifying Interaction Effects and Selection of Regression Models
The three different kinds of interactions in a spatial econometric model are divided into endogenous interaction effects between the explained variable (Y), exogenous interaction effects between the explanatory variables (X), and interaction effects between the error terms (ε). The spatial lag model (SLM) contains endogenous interaction effects, and its motive is that the urban land of a spatial unit depends on the urban land of adjacent spatial units. The mean value of urban land of adjacent spatial units, that is, the spatial lag value of the dependent variable, was incorporated into the regression as an explanatory variable. The spatial error model (SEM) contains the interaction effects between the error terms, and its motive is that the determinants of the explained variable omitted in the model are spatially related. Some unobservable factors jointly affect adjacent space units, but these factors are not included in the model, and the error term contains these common influencing factors, thus resulting in the spatial autocorrelation of the error term. The spatial Durbin model (SDM) includes endogenous and exogenous interaction effects, and its motive is that the urban land of a spatial unit depends on the urban land of adjacent spatial units, and also depends on the explanatory variables of adjacent spatial units. In other words, the spatial lag value of the explained variable and the spatial lag value of the explanatory variables are included in the regression as explanatory variables.
In the empirical study of the influencing factors of urban land, it is an important purpose to investigate the interaction effects. Starting from a standard linear regression model (LRM), the spatial lag model (SLM), spatial error model (SEM), and spatial Durbin model (SDM) were established. The interaction effects were identified by the significance of the coefficient ρ, λ, θ and the lag coefficients of the explanatory variables. In order to eliminate the influence of different dimensions of variables and alleviate the influence of heteroscedasticity on the estimated parameters of the model, the natural logarithms of the explained variable and the explanatory variable were taken before regression. This paper adopts the spatial weight matrix based on ROOK adjacency to calculate the spatial econometric models, and the results are shown in Table 2. In addition, the spatial weight matrix based on inverse Euclidean distance was substituted into the three models for calculation, and the results were compared with Table 2. The results show that although the coefficients of the estimated results are different, their signs and significance do not change fundamentally. This also indicates that the results are robust and reliable. Note: The value in parentheses is the t statistic of the coefficient; ***, **, and * mean passing the test at the significance level of 1%, 5%, and 10% respectively. The explained variable is Urbanland.
Firstly, linear regression model (LRM) was performed, and it was found that the global Moran's I of regression residuals was 0.163, and the Z score was 4.01. The Z score was greater than 2.58 and significant at the level of 0.01, indicating that the regression residual had significant spatial autocorrelation. This does not conform to the basic assumption that linear regression model (LRM) residual sequences are not autocorrelated, so the spatial regression model should be adopted.
WY is included in both SLM and SDM models. WY denotes the endogenous interaction effects among the explained variable. If the coefficient ρ of WY is not significantly 0, it indicates that there is an endogenous interaction effect among the explained variable. In this empirical study, the coefficient ρ was significant at the 0.01 level both in the SLM and SDM models. It indicates that the explained variable urban land scale has an endogenous interaction effect, and the local urban land scale is significantly affected by neighboring districts and counties. In SLM and SDM models, the values of ρ are 0.245 and 0.278, respectively. On the premise that other variables remain unchanged, an increase of 1% in the urban land scale of adjacent districts and counties will result in an increase of more than 0.245% in the local urban land scale. The scale of local urban land is significantly affected by neighboring districts and counties.
Wµ is included in the SEM model. Wµ denotes the interaction effects among the disturbance term of the different counties. If the coefficient λ of Wµ is not significantly 0, it indicates that there is an interaction effect among the disturbance term of the different counties. In this empirical study, the coefficient λ was significant at the 0.01 level in the SEM model. This means that the disturbance term satisfies E(µ i µ j ) 0, and the disturbance term of each district or county can be expressed as a function of the error term of adjacent districts or counties.
WX is included in the SDM model. WX denotes the exogenous interaction effects among the explanatory variables. If the coefficient θ of WX is not significantly 0, it indicates that there is an explanatory interaction effect among the explanatory variables. In this empirical study, the lag coefficients of the explanatory variables POP (total population) and Wages were significant at the 0.05 and 0.01 levels, respectively. This indicates that there is an exogenous interaction effect among the explanatory variables.
When selecting the model, we first compared LM (lag) and LM (error), and found that both the spatial lag model and the spatial error model were significant at the level of 0.01. Secondly, a comparison between the Robust LM (lag) and the Robust LM (error) showed that the Robust LM (lag) was significant at the level of 0.1, while the Robust LM (error) was not, so the SLM was suitable. However, an SLM only considered the spatial effects of the explained variable. The explanatory variables were also relevant. This circumstance was included in an SDM. The SDM model included both endogenous and exogenous interaction effects. In addition, compared with the R 2 and log-likelihood of the four models in Table 2, SDM model has the best fitting effect. Therefore, the SDM model was selected for subsequent analysis in this paper.

Decomposition Analysis of Spatial Effects
LeSage and Pace (2009) found that the point estimation of the coefficient may lead to erroneous conclusions [67]. They observed that the use of partial derivatives was more effective in testing the existence of spatial effects. By considering the feedback loop effects, the spatial effects are decomposed into direct effect, spatial spillover effect, and total effect [68]. The formula is as follows: After sorting out the formula, it can be obtained as follows: where S r (W) = V(W)(I n β + Wθ), V(W) = (I n − ρW) −1 . I n is the n order identity matrix, ln is the n × 1 order matrix whose elements are all 1, r = 1, 2, . . . , k, is the number of explanatory variables. Converted into the form of a matrix, we can get: 11 S r (W) 12 · · · S r (W) 1n S r (W) 21 S r (W) 22 · · · S r (W) 2n . . . . . . . . . . . .
S r (W) ij and S r (W) ii are the elements of the ith row and the i (or j) column in the matrix S r (W), respectively. S r (W) ij is the magnitude of the change in the explained variable of the ith space unit caused by the change in the explanatory variable X jr of the jth space unit. S r (W) ii is the magnitude of the change in the explained variable of the ith space unit caused by the change in the X ir of the ith space unit.
Since direct effect and spatial spillover effects are different for different spatial units, if there are N spatial units and K explanatory variables, K different matrices of N × N direct effects and spatial spillover effects will be obtained. Reporting these results succinctly is difficult. LeSage and Pace (2009) suggested that the direct effect should be measured by the mean value of the diagonal elements of the matrix, and the spatial spillover effect should be measured by the row or column mean of the off-diagonal elements of the matrix. The total effect is measured by the mean of the sum of all elements of the matrix. The spatial spillover effect can be obtained by subtracting the direct effect from the total effect.
Spatial spillover effects are composed of different coefficient estimates, and the degree of dispersion of spatial spillover effects depends on the degree of dispersion of all coefficient estimates involved. Therefore, it is impossible to judge whether the spatial spillover effect is significant based on the coefficient estimate and its standard error or t value. In order to judge the statistical significance of direct effects and spatial spillover effects, LeSage and Pace (2009) suggested using the variance-covariance matrix obtained from maximum likelihood estimation to simulate the distribution of direct effects and spatial spillover effects. The combination of D parameters can be obtained by sampling the variance-covariance matrix. The spatial spillover effect of a particular explanatory variable can be approximated by calculating the mean value of these D samples, and its significance level (t value) can be obtained by dividing its mean value by its corresponding standard deviation. µ kd represents the spatial spillover effect of the kth explanatory variable of the sampling d, then the overall spatial spillover effect of all samples and its corresponding t value will be: where µ kd is the spatial spillover effect of the kth explanatory variable, and t is its corresponding t value. According to the t value of the spatial spillover effect, we can test whether the kth explanatory variable has a spatial spillover effect.
Distinguishing the direct, spatial spillover effect, and total effect of socioeconomic variables helps to derive more accurate and unbiased results [68]. The empirical results of SDM decomposition analysis are shown in Table 3.  (1) Comparative Analysis of the Same Explanatory Variable The significance and sign of the direct effect and spatial spillover effect of the GDP in the secondary/tertiary sector have changed. The direct effect is significant while the spatial spillover effect is not. The direct effect sign is positive and the spatial spillover effect sign is negative. The elasticity coefficient of the direct effect of the GDP in the secondary/tertiary sector is positive, and it passed the test at the significance level of 0.01, indicating that the local urban land scale is significantly affected by the local GDP in the secondary/tertiary sector. Under the condition that other variables remain unchanged, the GDP in the secondary/tertiary sector increased by 1%, which will promote the local urban land scale to increase by 0.432%. The primary reason for the GDP in the secondary/tertiary sector in districts and counties is the increase in enterprise benefits. Therefore, enterprises may invest more capital to expand reproduction, and local governments provide more industrial land to support sustainable economic development, thus promoting the increase in urban land scale. The GDP in the secondary/tertiary sector in districts and counties indicates that the increase in the output value of the service industry inevitably requires land for public management and public services and commercial finance, which leads to an increase in urban land. When the GDP in the secondary/tertiary sector of neighboring districts and counties is large, it will attract local manpower, material and financial resources to neighboring districts and counties, and urban land is the carrier of economic activities. Therefore, the elasticity coefficient of the spatial spillover effect of the GDP in the secondary/tertiary sector is negative. This spatial spillover effect is not significant. This means that the growth of GDP in the secondary/tertiary sector of neighboring districts and counties has an inhibitory effect on the scale of local urban land, but the spatial spillover effect of this inhibitory effect is not obvious. The total effect of the GDP in the secondary/tertiary sector is positive but not significant. It shows that the economic development level of the BTHCDR does not have a prominent impact on the increase in urban land scale. With the upgrading of industries, the development of land-saving industries has become a trend.
The direct effect of the total population and spatial spillover effect were significant, and the signs were positive. The elasticity coefficient of the direct effect of the total population is positive, and it passed the test at the significance level of 0.01, indicating that the scale of local urban land is significantly affected by the local population. When other variables remain unchanged, the increase of the local total population by 1% will promote the increase in the local urban land scale by 0.287%. The larger the total population of districts and counties, the more residential land, green land, and public facilities land are needed to support it, which is bound to lead to an increase in urban land scale. The elasticity coefficient of the spatial spillover effect of the total population is positive, and it has passed the test at the significance level of 0.05, indicating that the spatial spillover effect of the total population is obvious.
When other variables remain unchanged, an increase of 1% in the total population of neighboring counties will promote an increase of 0.7% in the scale of local urban land. The total effect of the total population is positive and has passed the test at the significance level of 0.01. It shows that the impact of population growth in the BTHCDR on the increase of urban land is still very prominent, and controlling per capita construction land is still an important way to curb urban land growth in the future.
The significance and sign of the direct effect of urbanization rate and the spatial spillover effect have both changed. The direct effect is significant, while the spatial spillover effect is not. The direct effect sign is positive, while the spatial spillover effect sign is negative. The elasticity coefficient of the direct effect of urbanization rate is positive, and it has passed the test at the significance level of 0.01, indicating that the local urban land scale is significantly affected by the local urbanization rate. When other variables remain unchanged, the increase of local urbanization rate by 1% will promote the increase of local urban land scale by 0.021%. The improvement of urbanization rate in districts and counties means that the population is shifting from rural areas to cities. The new urban residents need more land for urban public services, such as education, medical care, and social security, and more land for infrastructure, such as roads, a centralized water supply, power supply, and heating. The expanded urban land scale provides a basic guarantee for people-oriented urbanization. The elasticity coefficient of the spatial spillover effect of the urbanization rate is negative, but it is not significant. That is to say, the growth of urbanization rate of neighboring counties has an inhibitory effect on the scale of local urban land, but the spatial spillover effect of this inhibitory effect is not obvious. The total effect of urbanization rate is positive, but only passes the test at the significance level of 0.1. The increase in urbanization rate and the advancement of urban system construction have led to an increase in urban land. The implementation of the new urbanization has promoted the integration of urban and rural areas, and the convenient facilities in rural areas have also been improved. The more livable rural areas have reduced the desire of some people to migrate to cities, making the total effect of the urbanization rate less significant.
The direct effect and spatial spillover effect of fixed-asset investment are not significant, and the sign is negative. Fixed-asset investment in local and adjacent districts and counties has an inhibitory effect on the scale of local urban land, but this inhibitory effect is not obvious. Fixed-asset investment occupies a huge amount of capital and affects the amount of capital needed for urban expansion and new areas. Therefore, fixed-asset investment has a certain inhibitory effect on the scale of urban land. The total effect of fixed-asset investment is negative and insignificant.
The direct effect of the average wages is not significant, while the spatial spillover effect is significant. The elasticity coefficient of the direct effect of the average wages is positive, but not significant. The increase in local average wages has a driving effect on the scale of local urban land, but this driving effect is not obvious. The spatial spillover effect of average wages is positive, and it passed the test at the significance level of 0.01, indicating that the local urban land scale is significantly affected by the average wages in neighboring counties. When other variables remain unchanged, an increase of 1% in the average wages of neighboring counties will promote an increase of 1.717% in the size of local urban land. At the district and county levels, the spatial spillover effect of average wages is obvious. A likely reason for this is that districts and counties with higher average wages, such as districts in Beijing, Tianjin, and prefecture-level cities in Hebei, still have limited purchasing power relative to local housing prices. They chose neighboring counties and districts with lower housing prices, such as Zhuozhou and Yanjiao in Hebei, and even Chongli, where the winter Olympics will be held. Demand from neighboring counties has helped boost land use in local cities and towns. The total effect of average wages is positive and has passed the test at the significance level of 0.01. Income is the source of all consumption. The increase in average wages in the BTHCDR and relatively complete social security have increased the purchasing power of urban residents. The increase in domestic demand has enabled the healthy operation of the entire economy, thereby increasing the demand for land in all walks of life.
The direct effect and spatial spillover effect of highway mileage density are not significant, and the signs are positive. The increase of highway mileage density in local and adjacent districts and counties has a pull effect on the scale of local urban land, but this pull effect is not obvious. With the increase of highway mileage density in each district and county, the accessibility between districts and counties can be improved, the speed of goods circulation can be accelerated, and the communication between districts and counties can be enhanced. The total effect of highway mileage density is positive, but only passes the test at the significance level of 0.1. Transportation integration is a pioneering area for the coordinated development of Beijing-Tianjin-Hebei, and road network construction will surely reshape the urban land layout and spatial structure.
The direct effect and spatial spillover effect of actual utilization of foreign capital are not significant, and the symbols are positive. Foreign investment has absorbed labor force employment, promoted labor force transfer from rural areas to urban areas, and raised the urbanization rate. In addition, foreign investment is mostly concentrated in the manufacturing industry, which improves the production technology, management concept and labor productivity of enterprises, promotes the development of the secondary and tertiary industries, and is an important factor driving economic growth. These two reasons indirectly promote the expansion of urban land, but not significantly. The total effect of actual utilization of foreign capital is positive and insignificant.
The direct and spatial spillover effects of topographic relief are not significant; the former is negative, and the latter is positive. The smaller the local topographic relief, the flatter the terrain, and the less constrained the increase in urban land. Smaller topographic relief is more conducive to the construction of "seven connections and one leveling", which is conducive to the increase in urban land. If the topographic relief of neighboring counties is relatively large, the population and industry may gather locally, which will promote the increase in local urban land. The total effect of topographic relief is negative and insignificant.
(2) Comparative Analysis of Contribution of Different Explanatory Variables In order to horizontally consider the difference of the influence of different explanatory variables on the explained variables, this paper aimed to clarify the contribution degree of different explanatory variables. First, the explanatory variables were graded according to their significance level, and then the elasticity coefficients of each explanatory variable were compared under the same significance level. This method combining the significance level and elasticity coefficient of explanatory variables could compare the contribution degrees of different explanatory variables to the explained variables in detail.
The order of the direct-effect elasticity coefficient is: the GDP in the secondary/tertiary sector > the total population > urbanization rate. The symbols of the elasticity coefficient are all positive. It can be seen that the increase in local GDP, total population, and urbanization rate will directly lead to an increase in local urban land scale. The fixed-asset investment, the average wages, the density of highway mileage, the actual utilization of foreign capital, and the topographic relief have no significant direct effect on the scale of local urban land. In terms of explanatory variable symbols, the elasticity coefficient of the average wages, the highway mileage density, and the actual utilization of foreign capital are positive, indicating that these three factors have a positive effect on the local urban land scale. The elasticity coefficient of fixed-asset investment and topographic relief is negative, which indicates that an increase in the local fixed-asset investment and topographic relief will negatively impact the local urban land scale.
The order of elasticity coefficient of spatial spillover effect is: the average wages > the total population > highway mileage density. The symbols of the elasticity coefficient are all positive. An increase in the average wages, the total population, and the density of highway mileage will promote an increase in the local urban land scale. The GDP in the secondary/tertiary sector, urbanization rate, fixed-asset investment, actual utilization of foreign capital, and topographic relief are not significant. In terms of the symbol of explanatory variable elasticity coefficient, the actual utilization of foreign capital and the elasticity coefficient of topographic relief degree are positive, indicating that these two factors of neighboring districts and counties have a positive effect on the local urban land scale.
The elasticity coefficients of the GDP in the secondary/tertiary sector, urbanization rate, and fixed-asset investment are negative, indicating that an increase in the GDP in the secondary/tertiary sector, urbanization rate, and fixed-asset investment of neighboring districts and counties will negatively inhibit the increase in the land scale of the local urban land scale.

Comparison with Traditional Cross-Sectional Regression
This paper applies the method of spatial econometrics to the study of urban land at the county level in the BTHCDR, making up for the lack of correlation between spatial units in previous studies. The spatial correlation test results of traditional regression residuals of cross-sectional data show that residuals have spatial autocorrelation. There is an obvious spatial dependence on the scale of urban land between the districts and counties in the BTHCDR. The urban land in each county shows not a completely random state, but a spatial relation between similar values. The lag term or spatial error term of the spatial endogenous variables of SLM, SEM, and SDM reached 0.245, 0.402, and 0.278, respectively, and were significant at the 0.01 level, which proves that there is an obvious spatial correlation between the growth of urban land among counties. Specifically, this spatial correlation is manifested as positive external spillover between districts and counties. That is, the growth of urban land in districts and counties will have an effect on the growth of urban land in neighboring districts and counties through space overflow. It is suggested that, in the previous theoretical and empirical studies of urban land, the traditional regression analysis that ignores spatial dependence is deficient in theory and inconsistent with reality. As the coordinated development of the BTHCDR has become a national strategy, it is necessary to break down administrative barriers and take the BTHCDR as a whole to optimize the distribution pattern and spatial structure of urban land. Combined with the characteristics of urban economic development, from the perspective of urban land supply and support for economic development, the urban land supply linkage between regions should be innovated. Regional coordinated development planning is formulated to provide a scientific basis for the integration of Beijing, Tianjin, and Hebei, and fully realize the spilt economic benefits of Beijing and Tianjin [46,69]. At the same time, we should not only focus on the growth of the absolute size of the city, but also look at the improvement in the quality of urban land. Urban land optimization should limit the increment and tap the stock potential as the main line. The new urban land supply is mainly non-capital function transfer and Beijing-Tianjin industry transfer. Stock tapping aims at increasing per capita output and encourages the conversion of urban stock land into ecological land to realize green and sustainable development.

Comparison with Other Studies
In order to compare these results with scholars' conclusions, the positive or negative effects of explanatory variables on urban land were firstly analyzed without distinguishing between direct effects and spatial spillover effects. GDP in the secondary/tertiary sector has a positive impact on urban land scale, which is consistent with previous research results [15,21,28,36,42]. GDP is an important indicator of the governance ability of local governments, which can obtain huge benefits from the conversion of agricultural land to urban land [21]. Total population has a positive impact on urban land scale, which is consistent with previous research results [15,21,28,37,42]. People are the subject of all production and life activities. Intuitively speaking, more people take up more space. The urbanization rate has a positive impact on urban land scale, which is consistent with previous research results [15,42]. Increased urbanization means the migration of people from rural areas to cities, which in turn requires more land to accommodate these growing populations. Average wages or income have a positive impact on the scale of urban land, which is consistent with previous research results [36,37,[70][71][72]. High-income areas will need to build more houses, shopping centers and public facilities, which will use more land. Actual utilization of foreign capital has a positive impact on urban land scale, which is consistent with previous research results [15,42]. The establishment of service or high-tech companies requires land [73]. The high density of highway mileage makes it easy to transport the goods of other areas and counties. The conclusion that highway mileage density has a positive effect on urban land scale is consistent with previous studies [28,42]. Topographic relief has a negative impact on urban land scale, revealing that urban expansion is restricted by natural conditions and a flat terrain is more conducive to urban land expansion, which is consistent with previous findings [21,28,42,70]. The impact of fixed-asset investment on urban land scale is negative, which is inconsistent with previous research results [42]. A possible reason for this is that a considerable proportion of fixed-asset investment is used for the renewal and reconstruction of urban built-up areas, which occupies a huge amount of capital and affects the amount of capital needed for urban expansion and new areas, thus restraining urban expansion. For districts and counties with large investment volume of fixed assets, such as Beijing, Tianjin, and Beijing the investment amount of fixed assets per urban land is relatively large, which indicates that urban land is more economical and intensive.
Different from previous studies, this paper compares the positive and negative changes of the direct effect and spatial spillover effect coefficients of each explanatory variable, which can provide a basis for the formulation of coordinated development policies in Beijing, Tianjin, and Hebei. The total population is the new force of all production activities and the source of domestic demand. Average salary determines people's quality of life. The density of highway mileage reflects the convenience of the connection between districts and counties. The actual use of foreign capital does not crowd out county fiscal revenue and increase investment. An increase in these variables will promote the increase of urban land scale in both local and neighboring counties. The influence of fixed-assets investment on the scale of urban land in both local and adjacent counties is negative. It shows that, in the new era of the coordinated development of Beijing, Tianjin, and Hebei, urban development should focus on its own connotation construction, taking saving land as the starting point, and promote urban sustainable development. The growth of local GDP and the urbanization rate has a positive influence on the land scale of local cities and towns, and a negative influence on the land scale of neighboring districts and counties. It shows that the agglomeration effect of economic activities and urbanization on urban land scale is greater than the diffusion effect. When the growth rate of the local economy and urbanization exceeds that of surrounding counties, the growth of urban land in surrounding counties will slow down. Therefore, against the background of the coordinated development of the BTHCDR and new-type urbanization, it is necessary to consider the relationship between county economic growth and urban land growth in the process of adjusting urban layout and optimizing urban scale and structure in the BTHCDR, so as to build an urban-rural spatial pattern coordinated with economic development.
In terms of direct effect contribution, the economy, total population and urbanization rate are the main influencing factors on the urban land scale in the research area, which is similar to previous research results [15,21,28,36,42]. The difference is that the three variables have different degrees of contribution. This may be due to differences in region, scale, or period. In the BTHCDR, the contribution of GDP in the secondary/tertiary sector is significantly higher than other factors at the county level, which means that the economy is the most important influencing factor, consistent with previous studies [42]. The secondary industry, that is, the factory area of the processing manufacturing industry, is mostly distributed in cities and towns, and large-scale commercial and financial industries and other tertiary industries are mostly distributed in cities and towns. Urban land is the carrier of economic activities and the place where the urban population survives. The contribution of total population to urban land scale is second only to the economy. Population is the inexhaustible driving force of all economic activities and the internal cause of the increase in urban land scale. In the BTHCDR, the urbanization rate contributes little to urban land scale at the county level. This indicates that the urbanization rate of most counties in Hebei is at a low level, and their influence on urban land scale is lower than that of the economy and population. In terms of the contribution of the spatial spillover effect, the most important finding of this paper, is that the spatial spillover effect of average wages is the most prominent influencer. Whether in economically developed districts and counties or economically backward districts and counties, the local average monthly salary is significantly lower than the local housing price per square meter, which causes the residents in counties with a higher average salary to buy houses in counties with lower housing prices. Demand from neighboring counties has spurred local developers to buy new land and build residential complexes. The urban expansion caused by the increase in residential land is the main contribution to the increase in urban land scale in the BTHCDR.

Conclusions
At the county level, three spatial econometric models are used to investigate the endogenous interaction effects, exogenous interaction effects, and the interaction effects between error terms. Based on the perspectives of direct effect, indirect effect, and total effect, this paper innovatively analyzed the contribution degree of socioeconomic factors to urban land scale in the BTHCDR. Three conclusions were summarized as follows: (1) Urban land has an endogenous interaction effect. The local urban land is significantly affected by neighboring districts and counties. An increase of 1% in the urban land scale of neighboring districts and counties will lead to an increase of more than 0.245% in the local urban land scale. (2) The local population, average wages, highway mileage density, and actual utilization of foreign capital have a positive effect on the scale of urban land in local and neighboring districts and counties. Local fixed-asset investment has a negative effect on the urban land scale of the local and adjacent districts and counties. Local GDP in the secondary/tertiary sector and urbanization rate have positive effects on local urban land scale, but negative effects on adjacent districts and counties. Local topographic relief has a negative effect on the scale of local urban land, but a positive effect on the scale of adjacent districts and counties. (3) The order of direct effect contribution degree of different explanatory variables is as follows: GDP in the secondary/tertiary sector > total population > urbanization rate. The signs of elastic coefficients are positive and significant. The direct effects of fixed-asset investment, average wages, highway mileage density, actual foreign capital, and topographic relief on urban land scale are not significant. The contribution degree of the spatial spillover effect of different explanatory variables is ranked as follows: average wages > total population > highway mileage density, and the signs of elasticity coefficient are positive and significant. The spatial spillover effect of the urbanization rate, GDP in the secondary/tertiary sector, fixed-asset investment, actual utilization of foreign capital, and topographic relief are not significant.
The theoretical significance of the research results of this paper is to provide a theoretical basis for the optimization of urban land distribution in the BTHCDR. In terms of land-use practice, the state implements strict policy of controlling the total amount of new construction land, and the indexes of land needed for urban construction are insufficient. So, how to solve the shortage of land-use indicators? It is undoubtedly an ideal way to give full consideration to the influence factors of urban land scale, rationally arrange the use timing sequence of new urban construction land in the cities of Beijing, Tianjin, and Hebei, and accurately place the construction land index in the first place that needs construction. Through the macro-control policy on the supply of new urban land, we can guide the orderly flow of resource factors, promote the development of complementary and mutually promoting industries, and then realize the coordinated development. Practical operability guidance is the most important original intention of this paper. However, the paper also has some limitations. It is difficult to consider all the influencing factors, and only eight influencing factors were selected for the current study. The contribution degree of influencing factors to urban land varies with time. Future studies may attempt to add other influencing factors as well as to take long-term data for research. The variables used in this study, such as urban land scale, GDP, and actual utilization of foreign capital, were absolute values rather than per capita. Compared with economically backward counties in Hebei province, the absolute value of GDP of municipal districts in big cities is relatively high. Whether this absolute value index will have corresponding impact remains to be further studied by scholars using the per capita value. The authors look forward to use the results to further analyze the forces at work in the urbanization processes in this region of China and the spatial planning conclusions that would be drawn.