Temporal and Spatial Variability of Carbon Emission Intensity of Urban Residential Buildings: Testing the Effect of Economics and Geographic Location in China

Highlights

A new method to calculate the URBCE is proposed for explore the influencing factors of carbon emission.

The spatial econometric model solved the bias estimation caused by spatial correlation.

The effect of regional factors on carbon emission is determined by an econometric model.

Regional emission reduction policies are more scientific and effective than national ones.

Abstract

The role of urban residential buildings (URBs) in the carbon reduction goal of China is becoming increasingly important because of the rising energy consumption and carbon emission of such buildings in the region. Considering the increasing spatial interaction of the carbon emission of URBs (URBCE) in the region, this study investigates the influence of climate and economic factors on the URBCE in North and South China. First, the URBCE is calculated by using a decomposition energy balance table based on the carbon emission coefficient of electric and thermal power, thereby improving the estimation of the basic data of URBCE. Second, the influence of economic and climatic factors on the URBCE intensity in 30 provinces of China is explored by using a spatial econometric model. Results show that the URBCE intensity in China had a spatial autocorrelation from 2000 to 2016. Climatic and economic factors have great differences in the degree and direction of influencing the URBCE intensity in the country. Formulating emission reduction policies for climate or economic zones is more scientific and effective than developing national policies. Among these factors, urbanization rate, climate, and GDP per capita have a significant positive impact on the URBCE intensity in the region, whereas other factors have varying degrees of negative impact. In addition, climate, consumption level, and building area have significant spatial spillover effects on URBCE intensity, whereas other factors do not pass the significance test. Relevant conclusions should be given special attention by policymakers.

1. Introduction

An International Energy Agency survey (IEA, 2019) found that the residential building sector is responsible for more than 20% of the final energy consumption worldwide. China ranks first in residential energy consumption in the world [1]. In the past decade, the amount of carbon dioxide released by China’s residential building sector has grown rapidly, with a 6.57% increase annually, and the carbon emission of urban residential buildings (URBCE) exceeded 1.2 billion tons in 2016 [2,3] (Figure 1). Coming from the second-largest energy consumer in China, the issue of URBCE cannot be ignored [4]. With the continuous advancement of the urbanization (UR) process in China [5], urban population and urban residential building (URB) areas are also growing, and the energy demand and URBCE will be increased further [6,7]. This sector will contribute to over 50% of energy savings needed to reach the goal to peak carbon emissions ahead of 2030 [8]. Therefore, effective mitigation of the URBCE in China is crucial to the realization of the Paris Agreement in 2030 [9,10].

Many scholars have made great contributions in the exploration of the main factors that affect URBCE. Nie et al. found that the improvement of the energy efficiency of household appliances could reduce the energy use of residential buildings [11]. Highly efficient energy microgrids, distributed energy [12,13], and zero-energy consumption building technology [14] could also be used to reduce energy waste. Many scholars have also adopted the decomposition of the factors that influence the energy consumption intensity of residential buildings, such as the LMDI model [15,16] and the STIRPAT model [17,18] to explore the effective ways of reducing carbon emission from residential buildings. However, the climatic environment and economic development levels among regions are quite different because of the vast territory of China, thereby resulting in large differences in the URBCE intensity among its provinces. Because of the regional climate factors determine the energy consumption model of residential buildings [19], such as heating in the north [20] and cooling in the south [21]. However, the aforementioned studies ignore the climatic and economic factors and spatial interaction effects on the carbon emissions of residential buildings [22,23]. To fill this research gap, Miao discussed the effects of climatic and economic factors on energy intensity and the URBCE based on the data on Chinese urban residents in 2013 and found that climatic factors are positively correlated with URB energy consumption [17]. However, the conclusion possesses uncertainty because of the short research time span. Nie et al. compensated for the shortcoming of the aforementioned research by supplementing the time series data and confirmed the important role of climatic factors in the energy consumption of residential buildings [11]. Although these studies strongly support regional carbon emission reduction policy making, the analysis of factors that affect the carbon emission of residential buildings is still not comprehensive enough to draw robust conclusions. For example, the spatial interaction effect of URBCE is often ignored in previous research and policy making because of the prominent “fixed” characteristics of URBs [24]. Consequently, similar research conclusions are flawed [25].

The research conclusion of a spatial econometric model achieved greater explanatory power than traditional models because it improves the performance of the panel data model and solves the biased estimation caused by spatial correlation [24]. Moreover, spatial econometric models are widely used to study the drivers of carbon emission [26,27,28]. For example, Zhu et al. investigated the effect of climatic factors on pollutant emissions based on a spatial econometric model of electricity power data and found that climate warming will reduce residential energy use fuel and increase electricity consumption [29]. However, the results in a given energy type, such as the sole use electric power data, are not thorough enough for generating robust conclusions. Zhao et al. analyzed the climatic and economic factors that affected the URBCE by using a spatial econometric model [30]. This model is mainly based on the fitting of nighttime light datasets and statistical carbon dioxide emissions to obtain new carbon emission data. The results show that the regions performed better than the nation in estimating the URBCE. Moreover, the spatial effects of climatic factors have an important impact on the carbon emissions of urban residents, thereby providing an important basis for the decision-making objectives of regional emission reduction [22]. However, the database contains biased estimates. The reason for this is that the nighttime lighting data used match the total carbon dioxide emissions of cities but do not reflect the residential carbon emissions (Appendix B). Therefore, the influencing factors of the spatial difference of URBCE intensity should be studied by a spatial econometric model, and regional emission reduction strategies must be formulated accordingly [29,30].

The first step, and a prerequisite of this research, is accounting for the URBCE data. In existing studies, the accounting methods of building energy consumption and carbon emissions are mainly divided into two methods: micro-energy consumption statistics and macro-energy consumption measurement. Micro-energy consumption statistics usually uses the terminal energy consumption of buildings to realize the statistics of building energy consumption and carbon emissions. Bastos et al., calculated the carbon emissions of three types of residential buildings using the life cycle method, in which the carbon emissions of electricity and heat contributed the most to the overall carbon emissions [31]. Scheuer et al. regarded a public building at the University of Michigan as the research object. Assuming that its life cycle energy consumption was estimated in 75 years, the results showed that HVAC and electricity accounted for 94.4% of the life cycle primary energy consumption [32]. Many macro measurement methods, such as the emission coefficient method (IPCC inventory method), input–output method [33,34], life cycle assessment method [35], and material balance algorithm [36], are available for building carbon emissions, each of which has its application scope and relative merits. In summary, most previous studies reveal a lack of a unified standard for building carbon emission calculations, and the final results of building carbon emission calculations vary greatly, thereby causing difficulties in providing a reference for the government’s plan for energy saving [37].

With this in mind, this study tries to fill these gaps, and aims to make the following contributions. First, this work constructs a regional carbon emission intensity calculation model that is superior to the whole country. It is used to calculate the carbon emission intensity of urban residential buildings in district heating and distribution heating zones, in which the carbon emission factors of electric power and thermal energy are distinguished and calculated separately. The calculation results provide scientific and accurate basis for energy conservation and emission reduction in the construction field, and propose corresponding emission reduction strategies according to regional characteristics [38]. According to the current situation of research on the spatiotemporal distribution of carbon emissions in the building sector, relatively few works exist on the spatiotemporal distribution characteristics of existing carbon emissions. Moreover, the overall situation and development trend of national or regional building carbon emissions are difficult to grasp because of the short research cycle and the fact that the research object of the existing literature is only the direct carbon emissions of the construction industry [37,39]. From the perspective of the spatiotemporal evolution mechanisms of the URBCE, a national accounting model of carbon emissions from existing buildings in provinces, cities, or regions should be established urgently at the macro level, and the spatiotemporal characteristics of carbon emissions (energy consumption types, carbon emission types, carbon emission intensity, and geographic space) must be revealed immediately [40]. Second, through a spatial econometric model, the influencing factors of the spatial difference of URBCE intensity are analyzed. The conclusions can help government departments in formulating local energy saving and emission reduction targets to promote the sustainable development of a low-carbon economy. The process of solving the problem is shown in Figure 2.

The rest of this paper is organized as follows: Section 2 provides the methodology of URBCE and spatial autocorrelation analysis. Section 3 presents the variable selection and data source. Section 4 presents the spatiotemporal evolution of carbon emission intensity. Section 5 presents the analysis on URBCE intensity and influencing factors. Section 6 presents the conclusions, policy implications, limitations, and future directions.

2. Method

2.1. Carbon Emissions and Intensities of URBs

In operation, URBCE mainly comes from various energy sources used by urban residents in their daily life, such as heating, air conditioning, lighting, cooking, and hot water and household appliances, which includes three types of fossil energy sources (e.g., coal, liquefied petroleum gas, and natural gas) and two types of secondary energy sources (e.g., electric and thermal power). Given the substantial difference in climate between the north and south regions, 15 regions in the north belong to the district heating area, whereas 15 regions in the south belong to the distributed heating area (Figure 3) [41]. According to the characteristics of energy consumption, the energy consumed in the northern region is divided into four categories, namely, liquefied petroleum gas, natural gas, electricity, and heat [42]. Similarly, the energy consumed in the southern region is divided into four categories, namely, coal, liquefied petroleum gas, natural gas, and electricity.

2.1.1. Carbon Emission in URBs

The calculation method of this study is mainly determined by the carbon emission factors of various types of energy and energy consumption of URBs. Among the carbon emission factors of various types of energy, the carbon emission coefficient of fossil energy can be obtained on the basis of calculation methods or direct data published by relevant institutions at home and abroad. However, the power carbon emission coefficient only has a few years of public data for reference, and nearly no public data on the thermal carbon emission coefficient are available for inquiry. This situation causes most researchers to only consider fossil energy emissions when calculating carbon emissions. The present study fully considers all the aforementioned carbon emission factors and constructs a carbon emission calculation model for URBs in the northern and southern regions.

$$E_N=C_n{f}_n+C_l{f}_l+C_e{f}_e+C_h{f}_h,$$

(1)

$$E_S=C_n{f}_n+C_l{f}_l+C_e{f}_e+C_c{f}_c,$$

(2)

where $E_N \left(kgCO2\right),E_S$ (kgCO₂) represent the carbon emissions of residential buildings in the northern and southern regions, respectively; $C_n\left(\mathrm{cubic}\mathrm{metre}-\mathrm{cu}.\mathrm{m}\right), f_n\left(kgCO2/\mathrm{cu}.\mathrm{m}\right)$ are the consumption and the carbon emission coefficient in different regions, respectively;$C_l\left(tons-t\right), f_l\left(kgCO2/t\right)$ are the consumption and the carbon emission coefficient of liquefied petroleum gas, respectively; $C_e\left(kwh\right), f_e\left(kgCO2/kwh\right)$ are the power consumption and the carbon emission coefficient of electric power, respectively; $C_h\left(MJ\right), f_h\left(kgCO2/MJ\right)$ are the thermal energy consumption and the carbon emission coefficients of the north region, respectively; $C_c\left(t\right), f_c\left(kgCO2/t\right)$ are the coal consumption and the carbon emission coefficients of the southern region, respectively.

2.1.2. Carbon Emission Coefficient of Electric Power

The consumption of electricity in the life of urban residents in China is gradually increasing because of the improvement of China’s economic development level and the continuous enhancement of people’s living standards. The carbon emission coefficient of electricity must be understood first to measure the carbon emissions generated by electricity consumption accurately. China’s regional grid baseline emission factors published by the National Development and Reform Commission were used in many early studies. These carbon emission factors are announced annually; they are measured by carbon emission reduction, which is mainly used for the development of CDM projects, thereby misleading many research institutions or scholars [43]. In recent years, relevant institutions have provided other types of power carbon emission factors or their measurement methods, which were issued by the “Guidelines for the Preparation of Provincial Greenhouse Gas Inventories” and the National Development and Reform Commission. According to the “Qin Mountains-Huai River” boundary, the Chinese government set district heating and distributed heating zones in 1957, named the North–South heating line. District heating zones include North China, including northeast and northwest regions of China. At the same time, due to the different efficiency of power generation and the different energy structure, the average carbon emission coefficient of each region is quite different (Figure 3). China’s power grid is divided into six regions (e.g., north, northeast, east, central, northwest, and south China) in accordance with the National Development and Reform Commission’s method of dividing the regional power grid to reflect the characteristics of power supply transfer among different regions and the average carbon emission coefficient of regional power grids. The carbon emission factors of the North-South zones heating mode and the regional power grid are the basis for calculating the carbon emission of buildings in various provinces and cities.

In order to obtain data easily, the data of fossil energy such as coal, oil, natural gas consumed in the power generation process and the power generation capacity of each province (including thermal power generation and renewable energy power generation) are all taken from the energy balance table of each region in the China Energy Statistics Yearbook. The average carbon emission coefficient of each regional power grid is presented as follows:

$$f_e=\left({{\displaystyle \sum }}_{i=1}^n{{\displaystyle \sum }}_{j=1}^m{C}_j{f}_j\right)/{{\displaystyle \sum }}_{i=1}^n\left(T_i+R_i\right),$$

(3)

where $f_e\left(\mathrm{kgCO}2/\mathrm{kwh}\right)$ is the average carbon emission coefficient of each regional grid; $C_j\left(t\right), f_j\left(kgCO2/t\right)$ are the energy consumption and the carbon emission coefficient of the j-th fossil in the i-th region, respectively; $T_i\left(kwh\right),R_i\left(kwh\right)$ are the amounts of thermal power generation and renewable energy generation, respectively. According to the formula (3), the national carbon emission coefficient in 2015 can be calculated to be 0.6131 kgCO₂/kwh. The national electricity carbon emission coefficient of the most recent year (2015) published by the National Development and Reform Commission is 0.6101 kg CO₂/kwh, which is close to 0.6131 kg CO₂/kwh. Therefore, the method is considered feasible. On the basis of the aforementioned formula, this study calculates the carbon emission coefficient of each region from 2000 to 2016 separately.

2.1.3. Thermal Carbon Emission Coefficient

Affected by climate differences, the heating demand in winter is strong in the cold areas of North China, especially in the northeast, where coldness lasts longer than five months and has a great impact on the URBCE. The thermal carbon emission coefficient must be determined to measure carbon emissions from thermal consumption accurately. However, the influence of energy efficiency and consumption ratio of coal-fired boilers, gas-fired boilers, and cogeneration equipment should be considered comprehensively. No accurate and reliable thermal carbon emission coefficient exists for reference because of the difficulty of obtaining various basic data. In this study, the carbon emission intensity method is used to calculate the thermal carbon emission coefficient in the heating areas in the northern region.

$$E{F}_h={{\displaystyle \sum }}_{i=1}^n\left[\left({{\displaystyle \sum }}_{j=1}^m\frac{e{e}_i\times e{t}_{ij}}{C_j}\times f_j\right)\times pro{p}_i\right]$$

(4)

$$f_h=E{F}_h/{{\displaystyle \sum }}_{i=1}^n\left(e{e}_i\times pro{p}_i\right)$$

(5)

where $E{F}_h\left(kgCO2/㎡\right)$ is the comprehensive carbon emission intensity, $e{e}_i$ is the energy consumption intensity of heating technical type (e.g., coal-fired boilers, gas-fired boilers, and cogeneration of heat and power), $e{t}_{ij}$ is the proportion of j-th energy consumed by the i-th heating technique, $C_j$ is the conversion coefficient of standard coal for the j-th energy, $f_j$ is the carbon emission coefficient of j-th energy, and $pro{p}_i$ is the area proportion heated by the i-th technique.

2.2. Spatial Autocorrelation Analysis

Given the prominent “fixed” characteristics of URB, previous studies often neglected the spatial interaction effect of URBCE and only analyzed the URBCE intensity in time series, thereby generating biased conclusions. Spatial autocorrelation analysis is the precondition of setting up the spatial panel measurement model correctly. This type of analysis can measure the concentration or dispersion of URBCE intensity and reveal the evolution of the spatiotemporal patterns of URBCE intensity in China. On the basis of research scope, spatial autocorrelation can be divided into global and local spatial autocorrelation. Global spatial autocorrelation reveals the degree of spatial autocorrelation of an entire region, whereas local spatial autocorrelation reflects the degree of correlation between each spatial unit and its neighbors on a certain attribute [41,44].

2.2.1. Global Spatial Autocorrelation

This work uses the global spatial autocorrelation index (Global Moran’s I [GI]) to reflect the spatial autocorrelation characteristics of the overall URBCE intensity in the 30 provinces. GIS mapping techniques (especially ArcGIS software) can help to identify spatial autocorrelation provinces visually, but not statistically [45].

$$GI=\frac{n{{\displaystyle \sum }}_{i=1}^n{{\displaystyle \sum }}_{j=1}^n{W}_{ij}\left(y_i-\overline{y}\right)\left(y_j-\overline{y}\right)}{{{\displaystyle \sum }}_{i=1}^n{{\displaystyle \sum }}_{j=1}^n{W}_{ij}{{\displaystyle \sum }}_{i=1}^n{\left(y_i-\overline{y}\right)}^2},$$

(6)

where n is the total number of regions (n = 30); $y_i, y_j$,y are the carbon emission intensities of regions i and j, respectively; $\overline{y}$ is the average carbon intensity (LCI) of the residential buildings in all regions; and $W_{ij}$ is the spatial weight matrix, which reflects the proximity of spatial elements. The first-order rook weight matrix is used. The range of GI is between −1 and 1. The closer GI is to 1, the higher the agglomeration degree of spatial units with similar attributes; the closer GI is to −1, the higher the agglomeration degree of spatial units with different attributes. When GI = 0, spatial units do not exhibit an agglomeration phenomenon on a certain attribute and are randomly distributed [44]. In this study, GI is tested by the Z value, and zero hypotheses are accepted or rejected according to the Z value on the basis of setting the significance level.

2.2.2. Local Spatial Autocorrelation

The GI scatterplot can reflect the degrees of association and correlation among local regions but not the significance of the degree of association in local regions [46]. Therefore, this research uses a combination of GI scatterplots and Local Indicators of Spatial Association (LISA) cluster plots to analyze the spatiotemporal pattern evolution of URBCE intensity in the 30 regions of China from 2000 to 2016. Quadrants 1 and 3 and 2 and 4 of the GI scatterplot represent the positive and negative spatial correlation between samples, respectively. The first, second, third, and fourth quadrants indicate that the area with high, low, low, and high observation values are surrounded by a high-value agglomeration area (HH), high-value area (LH), low-value area (LL), and low-value area (HL), respectively [41]. LISA agglomeration maps represent different spatial autocorrelation types by using various colors and can reflect the significance of the degree of correlation in local areas.

2.3. Spatial Panel Measurement Model

Traditional econometric model assumes the spatial independence of the sample unit. Thus, the estimation results cannot truly reflect the spatial correlation of the panel data, and the important information of the panel data is lost. Spatial econometrics is based on traditional econometrics to introduce the spatial interaction of panel data into a model by establishing a spatial weight matrix, thereby fully reflecting the “panel” characteristics of the data.

According to the different ways of reflecting spatial interaction effects, three main models of spatial econometrics, namely, spatial lag model (SLM), spatial error model (SEM), and spatial Durbin model (SDM), exist.

(1) SLM

If the spatial correlation of the panel data is derived from the selected factor variable, then it can be analyzed by introducing the spatial lag term of the dependent variable.

$$y_{it}=\delta {{\displaystyle \sum }}_{j=1}^n{w}_{ij}{y}_{jt}+\beta x_{it}+{\mu }_i+{\lambda }_t+{\epsilon }_{it},$$

(7)

where i represents the region (i = 1, 2, 3,..., N), t is the time (t = 1, 2, 3,..., T), δ is the spatial autoregressive coefficient, $y_{jt}\left(kgCO2/㎡\right)$ is the URBCE intensity in year t of region i, $x_{it}$ is the selected explanatory variable and indicates the observation value of the explanatory variable in year t of region i, ${\mu }_i$ is the spatial fixed effect, ${\lambda }_t$ is the time fixed effect of the model, and ${\epsilon }_{it}$ is a random error term.

(2) SEM

Certain variables closely related to the interpreted variables may be omitted when selecting the interpreted variables of the model. Thus, the spatial spillover effects of random errors may exist. SEM is expressed to deal with this issue as follows:

$$y_{it}=\beta x_{it}+{\mu }_i+{\lambda }_t+{\varphi }_{it},{\varphi }_{it}=\rho {{\displaystyle \sum }}_{j=1}^n{w}_{ij}{\varphi }_{jt}+{\epsilon }_{it},$$

(8)

where ${\varphi }_{it}$ represents the spatial error autocorrelation error term, and ρ is the spatial autocorrelation coefficient of the error term.

(3) SDM

If endogenous interactions, exogenous interactions, and autocorrelation of error terms exist simultaneously, then SLM and SEM are inapplicable. SDM needs to be set.

$$y_{it}=\delta {{\displaystyle \sum }}_{j=1}^n{w}_{ij}{y}_{jt}+\beta x_{it}+{{\displaystyle \sum }}_{j=1}^n{w}_{ij}{x}_{ijt}\gamma +{\mu }_i+{\lambda }_t+{\varphi }_{it},$$

(9)

$${\varphi }_{it}=\rho {{\displaystyle \sum }}_{j=1}^n{w}_{ij}{\varphi }_{jt}+{\epsilon }_{it},$$

(10)

where $\gamma$ is the spatial lag and explains the coefficient of the variable. In this study, the traditional mixed panel models, namely, SLM, SEM, and SDM, are used to conduct statistical tests. The optimal model is selected by comparing and analyzing the test results of each model [47,48]. The basic steps are presented as follows.

(1) The data are estimated and tested by the traditional hybrid panel model without considering the spatial interaction effect. Then, the residuals are checked for spatial autocorrelation. On this basis, the applicability of SLM and SEM to the problem is investigated, and the fixed effects that need to be included are evaluated. The least squares method (OLS) is used to estimate the traditional mixed panel model. The spatial autocorrelation of residuals is tested by using the Lagrange multiplier (LMlag, LMerror) and its robust form (R-LMlag, R-LMerror).

(2) According to the results of LMlag, Lmerror, R-LMlag, and R-LMerror, the applicability of SLM and SEM is evaluated. If the results of LMlag, Lmerror, R-LMlag, and R-LMerror statistics are significant, then SLM and SEM are applicable to the problems to be studied. However, SLM and SEM are non-nested models.

(3) Whether the nested SDM model can be simplified to SLM and SEM should be verified. The null hypothesis of which SDM can be simplified to SLM and SEM is $H_0:\gamma +\delta \beta =0$ and $H_0:\gamma =0$. The Wald test is used herein. If the result passes the 5% significance level, then SDM cannot be reduced to SLM or SEM; otherwise, it can be simplified to SLM or SEM.

3. Variable Selection and Data Source

3.1. Variable Selection

This study studies the influencing factors of URBCE intensity in 30 regions in China (Appendix A). Most of the research on the factors that affect the URBCE is concentrated on the carbon emissions of the entire country or a single area. Research on the spatial difference of carbon emission panel data for URB is limited. In comparison with the carbon emissions of the industry and the construction industry, those of the URBCE have distinct characteristics. Residential buildings are crucial for people’s daily life. Their main function is to provide a suitable living environment to meet the daily basic needs of people. Therefore, carbon from residential buildings are emitted to shape a suitable living environment and maintain the energy consumed in people’s daily life (Figure 4).

The factors that affect the carbon emissions of residential buildings include the following.

(1) Climate zone

The first part of URBCE is due to the consumption of energy by household heating or refrigeration [49]. China’s territory spans many climatic regions, and the heating and cooling methods of URBs are greatly affected by climate. In the cold areas of North China, URBs adopt centralized heating to secure the heating demand of residents. In South China, where the climate is warmer and the temperature in winter is generally above 0 °C, no district heating is adopted. According to the survey and statistics of the China Building Energy Conservation Association, from 2000 to 2016, heating energy consumption in the north heating area accounted for 70% to 85% of the total residential building energy consumption [11,29]; whereas the centralized heating energy consumption in the south accounted for 0%, and household energy consumption was mainly electric energy consumption in daily life. In this study, the heating degree days (HDD) is used to measure the difference in climate zones [11]. The main reason is that the matching degree between climate zoning and HDD in the different regions of China is relatively high, which can effectively measure climate division [50], so HDD is selected as the characteristic variable of climate zone.

(2) Economic development level

The level of economic development has a dual impact on the URBCE [17]. First, the higher the level of regional economic development, the higher the energy consumption and material demand of residents for heating and refrigeration, which will increase carbon emissions [51]. Second, regional economic development has improved, residents are interested in healthy and environmentally friendly lifestyles, and their awareness of environmental protection has been strengthened, thereby affecting their behavior and reducing carbon emissions. The per capita gross domestic product (PGDP) and urban residents’ consumption level (UCL) are used to represent the regional economic development level to truly reflect the impact of regional economic development level on the total carbon emissions and intensity of URB [17].

(3) Life consumption

Consumption in daily life is the main source of URBCE. The carbon emissions from the daily life of residents mainly come from the electric power consumption for household appliances and the use of gas for cooking and bathing. The impact on electricity consumption is mainly related to the number of household appliances and habits of use. Gas consumption of is closely related to the household population, household eating quantity, and bathing habits. The living habits of residents are difficult to quantify, and such data are difficult to collect. Thus, this work selects the household electrical appliances per household (HEA) and the energy consumption structure of URBs (electricity accounts for total energy consumption of household proportion [LES]) as the indicators of the impact of domestic consumption on the URBCE to characterize the impact of residents’ living consumption on the carbon emissions of residential buildings [52,53].

(4) Control variables

The control variables that affect the URBCE mainly include UR rate (UR), URB area (LA), and urban per capita living building area (PLA).

UR has a dual impact on the LCI of residential buildings. On the one hand, the increase in UR is accompanied by large-scale urban construction and rise in urban population [17,52]. Urban per capita living energy consumption exceeds 80% of per capita energy consumption in rural areas [50], thereby resulting in a sharp rise in total building energy consumption. On the other hand, areas with high URs tend to have increased levels of urban construction. People may increase investment in the R&D of environmental technologies, purchase intermediate equipment that is conducive to environmental protection, and promote the reduction of the energy consumption of residential buildings. Therefore, this study selects UR as the control variable of URBCE intensity [53].

An increase in LA will greatly raise the total URBCE (LC) and positively affect the LC of urban residential construction. However, in the case of a certain total floor area, an increase in PLA will decrease the energy consumption per unit of the building area and create a “dilution effect” similar to the concentration of brine or sugar water. This research selects LA and PLA as the control variables of URB LCI to reflect the impact of LA on carbon emission intensity fully. The selected influencing factors and quantitative indicators of URBCE intensity are shown in

Table 1. Influencing factors and quantitative indicators of URBCE intensity.

Category	Indicator	Code
Climate zone	Heating degree days	HDD
The economic level	Per capita GDP	PGDP
	Urban residents’ consumption level	UCL
Living consumption	The average number of household appliances per 100 households (air conditioning, TV, refrigerator, washing machine, shower water heater)	HEA
Control variable	Urban residential building energy consumption structure	UR
	(electricity accounts for the total energy ratio)	PLA
	Urbanization rate	LA

[54].

3.2. Data Source and Processing

The main research scope is China except the data of Xizang, Hong Kong, Macao, and Taiwan. By regarding the administrative areas of provinces and cities as the research object, this study divides the climate area (district heating area) and economic area (distributed heating) according to the different heating types and economic attributes. Data sources include government statistics and research institutions. See

Table 2. Description of the materials used for estimating CO₂ emissions.

Data	Data Description, Unit	Year	Source
Energy consumption	Fossil energy, Electricity, building heating	2000–2016	China Energy Statistical Yearbook, China City Statistical Yearbook, (http://data.stats.gov.cn) China Building Energy Consumption Report [2]
HDD	Heating degree days (HDD), denoting the sum of mean daily temperature below 18 °C in a year, °C·d	2000–2016	China Meteorological Administration (http://data.cma.cn)
PGDP	PGDP = (GDP/GDP index)/Permanent Resident Population, Chinese Yuan (CNY)	2001–2017	China Statistical Yearbook
UCL	Urban Resident Consumption Level (UCL) = UCL/UCL Index, CNY	2001–2017	China City Statistical Yearbook, China Statistical Yearbook
HEA	The average household electrical appliances ownership (HEA) per 100 households	2001–2017	China City Statistical Yearbook
UR, PLA, LA	Control variable	2001–2017	China City Statistical Yearbook

Note: Yuans have been converted to 2000 constant prices for all price-related variables.

for the specific data sources and descriptions.

4. Spatiotemporal Evolution of Carbon Emission Intensity

4.1. Temporal Evolution of LC and Its Intensity

The LC and LCI trends of URBs in China during 2000–2016 are shown in Figure 5.

In Figure 5, the total URBCE in China from 2000 to 2016 showed a continuous growth trend; that is, it increased from 313 million tons in 2000 to 721 million tons in 2016. From the average annual growth rate of the URBCE, the average growth rate from 2000 to 2016 was 5.35%. However, the annual average growth rate varied greatly at different stages and could be divided into three stages as a whole. During the “10th Five-year Plan (2001–2005 FYP)” period, the average annual growth rate was 10.41%, of which the growth rate in 2005 was the highest during the observation period (reaching 15.27%). In the “11th FYP (2006–2010)” period, the average annual growth rate slowed down from the “10th FYP” period to 4.36%. During the “12th FYP (2011–2015)” period, the average annual growth rate slowed further than that in the “11th FYP” period, decreasing to 1.44%, which was approximately 86% lower than the average annual growth rate during the “10th FYP” period. In particular, the total URBCE fell from 680 million tons in 2012 to 667 million tons in 2013 and became negative for the first time.

Figure 6 shows that the carbon emissions per unit area of URBs can be roughly divided into two stages. During 2000–2010, the carbon emissions per unit area of the URBs showed a small fluctuation, which was basically maintained at approximately 32.14 kg CO₂/m², of which the maximum value during the observation period in 2010 was 33.22 kg CO₂/m². In 2010–2016, the carbon emissions per unit area of the URBs showed a continuous downward trend. By 2016, this value fell to the minimum during the observation period at 25.89 kg CO₂/m², which was 22% lower than that in 2010. The annual growth rate of carbon emissions per unit area of the URBs was −0.89%; that is, the carbon emissions per unit area of the URBs were weakened entirely. In 2001 and 2003, the growth rates were the largest (approximately 5%). In 2013 and 2016, the decline rate was the largest (approximately 6.8%).

4.2. Regional Carbon Emissions and Intensity of URBs

The characteristics of carbon emissions and intensity of URBs in China’s provinces should be thoroughly understood to formulate differentiated energy saving and emission reduction strategies and thus effectively control the total amount and URBCE intensity in China. The carbon emissions and intensity of the URBs in the provinces in 2016 are shown in Figure 7.

Figure 7 indicates significant differences in the carbon emissions and intensity of URBs in China. First, the north heating areas generally have high carbon emissions; Shandong (85.05 million tons), Hebei (57.01 million tons), and Liaoning (51.46 million tons) are the top three emitters, and the bottom three are Hainan (1.57 million tons), Yunnan (4.85 million tons), and Ningxia (6.55 million tons). Data comparison reveals that although Shandong and Hebei are the top carbon emitters, the difference in their carbon emissions is 28.04 million tons. The difference in carbon emissions between Shandong (ranked first) and Hainan (ranked last) is 54 times, indicating a considerable difference among the provinces.

The regional distribution of carbon emissions per unit area of the URBs shows that the carbon emissions per unit area of the URBs in the provinces are greatly affected by climate and economy. The north heating areas generally had higher carbon emissions per unit area than the non-heating areas, especially in severely cold areas. In the non-heating areas, coastal areas with good economic development generally had higher carbon emissions per unit area than those in the central and western regions. The top three regions are Heilongjiang (72.19 kg CO₂/m²), Qinghai (67.12 kg CO₂/m²), and Jilin (64.74 kg CO₂/m²), and the last three are Yunnan (5.62 kg CO₂/m²), Jiangxi 7.21 kg (CO₂/m²), and Hunan 8.42 kg CO₂/m²). The carbon emissions per unit area in Heilongjiang Province are approximately 13 times that of Yunnan Province.

4.3. Spatial Evolution of Carbon Emission and Its Intensity

(1) Global spatial autocorrelation analysis

The GI of URBCE intensity in the 30 regions of China from 2000 to 2016 was calculated using the software ArcGIS, and its significance was tested by constructing a normal distribution by the random permutation method. The results are shown in Figure 8.

Figure 8 shows the remarkable spatial autocorrelation characteristics of URBCE intensity in the 30 regions from 2000 to 2016. The GI of URBCE intensity in China was always approximately 0.5 during the study period, and its positive statistics Z(I) passed the 1% significance test. The GI declined slightly from 2000 (0.58) to 2004 (0.51), indicating that the spatial concentration of LCI was decreasing during this period. In 2005–2009, the GI showed an upward trend; it was 0.57 in 2009, showing that the spatial concentration of LCI increased during this period. In 2009–2013, the GI reduced significantly, indicating that the spatial agglomeration characteristics of URB LCI decreased significantly in this stage. In 2013–2016, the GI was constant at 0.49, and the spatial autocorrelation of LCI in the URBs in this stage was stable. Overall, the spatial autocorrelation of carbon emission intensity of the residential buildings in the 30 urban areas in China presented a variation of local fluctuations and a small overall decline, but the overall spatial autocorrelation was strong [41].

(2) Local spatial autocorrelation analysis

At the same time, ArcGIS software can help to identify scatterplot of the GI of URBCE intensity visually (shown in Figure 9). For 17 years, that is, from 2000 to 2016, the spatial correlation of the carbon emission intensity of the residential buildings in different regions and towns did not change significantly. A total of 12 areas fell into the first quadrant (HH), and 14 areas fell into the third quadrant (LL). In addition, Sichuan was in the second quadrant, Shandong was in the fourth quadrant, Qinghai was in the vicinity of the axis of symmetry between the first and fourth quadrants, and Henan was in the vicinity of the axis of symmetry between the first and second quadrants. From 2000 to 2016, the spatial agglomeration characteristics of URBCE intensity in China remained stable; 40% of the areas belonged to HH, and 46% belonged to low-value agglomeration areas [41].

The first quadrant represents a region with high total energy consumption of URBs, and the surrounding neighborhoods have relatively high energy consumption. The second quadrant represents an area where the total energy consumption of URBs is low and surrounded by areas with high energy consumption in the surrounding cities. The third quadrant represents areas where the total energy consumption of URBs is low and the total energy consumption of urban buildings in surrounding areas is low. The fourth quadrant is a region surrounded by areas where the total energy consumption of URBs is relatively high and surrounded by the surrounding cities. Among them, the first and third quadrants show the aggregation effect of high-valued aggregation and low-valued aggregation, which is a reflection of positive spatial correlation; the second and fourth quadrants are spatial negative correlation (spatial outliers), which is expressed as the research object. There is spatial heterogeneity among different regions.

Therefore, in the quantitative proof of spatial correlation, the energy consumption of URBs in China is mainly spatially dependent, but there is also a small amount of spatial heterogeneity.

A GI scatterplot cannot reflect the significance of the regional correlation degree. Thus, this study further analyzed GI scatterplot with a LISA agglomeration map (Figure 10). LISA can further judge the correlation between the attribute values of the studied regional units and the attributes of neighboring regional units, and determine whether it belongs to spatial clustering or spatial isolation.

Figure 10 depicts that the 30 areas can be divided into two groups: high-value aggregation areas (HH) and low-value aggregation areas (LL). In the 17 years, Inner Mongolia, Jilin, Liaoning, Shanxi, Gansu, Ningxia, and Hebei (in Northwest, North and Northeast China) always belonged to HH. Except in 2003, Heilongjiang was an HH. In 2000, Beijing showed high-value aggregation, and no significant spatial correlation was observed for the rest of the year. Tianjin showed high-value aggregation in 2000 and 2009 but no significant spatial correlation in the other years. Some provinces in Southwest, South and Southeast China showed a spatial correlation of low-value agglomeration. URBCE intensity in Guizhou, Hunan, Guangdong, and Jiangxi was always characterized by low-value aggregation from 2000 to 2016. Except for 2012, a spatial correlation of low-value agglomeration was observed in Yunnan. Zhejiang had no significant spatial correlation before 2006 and began to show the spatial characteristics of low-value agglomeration after this year. The spatial correlation of Guangxi in 2009–2012 was insignificant, but the spatial characteristics showed low-value agglomeration in the rest of the year.

The GI scatterplot and LISA agglomeration plot shows that the spatial correlation of URBCE intensity in the 30 provinces from 2000 to 2016 was mainly HH and LL. No regional agglomeration characteristic of LH and HL was observed. The HH was mainly distributed in the Northwest, North, and Northeast China, whereas the LL was distributed in parts of Southwest, South, and Southeast China. From 2000 to 2015, the total number of regions that show spatial correlation was stable, accounting for 43% to 53% of the total number of regions in the scope of this study [41]. Other regions did not show a significant spatial correlation. The total number of HH ranged from 6 to 9, accounting for 20% to 30% of the total number of research areas, whereas the total number of LL ranged from 5 to 8, accounting for 17% to 27% of the total number of research areas. In addition, no transition occurred in the spatial correlation characteristics of each year.

5. Influencing Factors of URBCE Intensity

Analysis of the spatiotemporal evolution of URBCE intensity reveals that the URBCE intensity in the 30 regions had a noticeable spatial correlation. Therefore, this work selected a spatial measurement model to analyze the factors that influence URBCE intensity. According to the model selection strategy, the traditional mixed panel data without spatial interaction effects were estimated and tested for residuals. The results are detailed in

Table 3. Estimation results of the traditional mixed panel data model.

	No Fixed	Spatial Fixed	Time Fixed	Bidirectional Fixed
HDD	0.16984 ***	0.08208 *	0.15575 ***	0.01856
PGDP	0.034428 ***	0.005748	0.039358 ***	−0.00398
UR	0.389262 ***	−0.005292	0.314348 ***	−0.011398
PLA	−0.60038 ***	−0.574695 ***	−0.981907 ***	−0.58509 ***
UCL	0.004239	0.05117 ***	−0.01385	0.03659 **
HEA	−0.17904 **	0.735398 ***	−0.285278 ***	0.16474 *
LA	−0.06281 ***	−0.149296	−0.02699 *	−0.54928 ***
LES	−0.98947 ***	−0.39759 ***	−0.95077 ***	−0.431366 ***
R2	0.922	0.348	0.935	0.408
Corrected R2	0.921	0.339	0.934	0.400
LMlag	29.8875 ***	5.1126 **	2.4661	0.0024
R-LMlag	21.6022 ***	24.5803 ***	0.0013	13.4154 ***
LMerror	55.3915 ***	25.3481 ***	14.3463 ***	3.4558 *
R-LMerror	47.1062 ***	44.8158 ***	11.8815 ***	16.8688 ***

Note: The numbers in parentheses are the t statistics (* p < 0.10, ** p < 0.05, *** p < 0.01).

. LMerror and R-LMerror passed the 1% significance test. Therefore, the model residual has spatial autocorrelation, and SLM and SEM are superior to the traditional mixed panel data regression without spatial interaction effects. However, SLM and SEM are non-nested models. The applicability of SDM and whether SDM can further degenerate into SLM and SEM need to be assessed further. The results are shown in

Table 4. Spatial Durbin model (SDM) applicability test.

	No Fixed	Spatial Fixed	Time Fixed	Bidirectional Fixed
Wald spatial lag	167.54 ***	167.54 ***	167.54 ***	167.54 ***
LR test spatial lag	144.56 ***	144.56 ***	144.56 ***	144.56 ***
Wald spatial error	151.21***	151.21 ***	151.20 ***	151.21 ***
LR test spatial error	130.66 ***	130.66 ***	130.66 ***	130.66 ***

Note: The numbers in parentheses are the t statistics (*** p < 0.01).

.

Table 4 shows that the Wald statistics coefficients passed the 1% significance test. Therefore, SDM cannot be simplified to SLM or SEM. The Hausman test passed the 1% significance test; thus, this research selects a fixed effect SDM to analyze the spatial factors of URBCE intensity. The results in

Table 5. SDM estimates.

	Time Fixed		Spatial Fixed		Bidirectional Fixed
	Coefficient	T	Coefficient	T	Coefficient	T
HDD	0.189 ***	12.857	−0.024	−0.566	0.0055	0.125
PGDP	0.045 ***	5.669	0.018 ***	2.707	0.0003	0.045
UR	0.254 ***	4.259	−0.116	−0.928	0.0609	0.484
PLA	−0.909 ***	−8.19	−0.862 ***	−6.284	−0.621 ***	−3.915
UCL	−0.051 ***	−2.62	0.0059	0.322	0.016	0.838
HEA	−0.174	−2.161	0.169 *	1.718	0.073	0.725
LA	−0.027	−1.839	−0.09	−0.871	−0.544 ***	−4.323
LES	−0.87 ***	−29.48	−0.458 ***	−11.29	−0.471 ***	−11.69
W*HDD	0.112 ***	2.915	0.186 ***	2.74	0.231 **	2.48
W*PGDP	−0.002	−0.086	−0.065 ***	−3.22	−0.116 ***	−5.41
W*UR	0.204	1.474	0.479	1.551	1.189 ***	3.670
W*PLA	0.112	0.489	−0.053	−0.199	0.410	1.206
W*UCL	−0.331 ***	−6.679	0.0305	0.611	0.134 **	2.441
W*HEA	1.247 ***	8.299	0.202	1.454	−0.265	−1.26
W*LA	−0.244 ***	−6.44	0.352	1.471	−0.917 ***	−2.97
W*LES	0.192 **	2.325	0.46 ***	5.757	0.289 ***	3.294
R2	0.9519		0.9845		0.9856
Adj–R2	0.9513		0.5393		0.4971

Note: The numbers in parentheses are the t statistics (* p < 0.10, ** p < 0.05, *** p < 0.01).

show that the estimation results of the time-fixed effect SDM were better than those of the spatial-fixed effect and bidirectional-fixed effect models [55,56].

The estimation cannot measure the influence of the feedback effect effectively because of the introduction of the spatial lag and interpreted variables in SDM [56]. Therefore, this study evaluates the influence of the independent variables on the dependent variables further according to the direct and indirect effects of the independent variables. The direct and indirect effects of time-fixed SDM are shown in

Table 6. Direct and indirect effects of time-fixed SDM.

	Direct	T	Indirect	T
HDD	0.194 ***	13.101	0.156 ***	3.643
PGDP	0.045 ***	5.732	0.005	0.211
UR	0.264 ***	4.174	0.271	1.679
PLA	−0.906 ***	−7.806	−0.019	−0.073
UCL	−0.063 ***	−3.304	−0.382 ***	−6.556
HEA	−0.135	−1.64	1.380	7.769
LA	−0.036 **	−2.364	−0.281 ***	−6.008
LES	−0.867 ***	−29.16	0.078	1.052

Note: The numbers in parentheses are the t statistics. (** p < 0.05, *** p < 0.01.)

.

(1) The direct effects of UR, HDD, and PGDP were 0.264, 0.194, and 0.045, respectively, indicating that the increase in these variables by 1% will improve the URBCE intensity in the region by 0.26%, 0.19%, and 0.045%, respectively. UR had the greatest impact, whereas PGDP had the least impact on URBCE intensity [52]. However, the indirect effects of UR and PGDP did not pass the significance test; only the indirect effect of HDD passed (0.156). Therefore, UR and PGDP had no significant impact on URBCE intensity in adjacent areas, whereas HDD changes can affect the URBCE intensity in such areas.

The above results show the outstanding fixed characteristics of URBs; that is, the non-mobility of residential buildings leads to reduced “communication” between the carbon emissions of the residential buildings among regions. In other words, the level of UR and economic development in adjacent areas will not have a significant spillover effect on the carbon emissions of the residential buildings in the region. In contrast, the spatial lag term of the HDD passed the 1% significance test, indicating that climate change in adjacent areas will have an impact on the LCI of the residential buildings in the area [50]. Climate had an obvious spatial similarity. Climate change in adjacent areas will have a significant impact on the climate of a region, which will lead to changes in the LCI of residential buildings in the region [11,29].

(2) The direct effects of PLA, UCL, LA, and LES were negative and passed the 5% significance test. PLA and LES had the greatest impact on the URBCE intensity in this area. The PLA and the proportion of electricity consumption in the URBs increase by 1%, which will decrease the URBCE intensity in this region by 0.9% and 0.87%, respectively. However, the indirect effects of PLA and LES did not pass the significance test, indicating that PLA and LES had no remarkable spatial spillover effect. The direct and indirect effects of UCL (−0.063 and −0.036, respectively) and LA (−0.382 and −0.281, respectively) on URBCE intensity passed the 5% significance test. Hence, the consumption level of urban residents and the area of URB will not only reduce the URBCE intensity in this region but also cause the URBCE intensity in adjacent areas. The consumption level of urban residents and the LA represent the regional development heat. According to the large current urban agglomerations in China, urban construction and economic development have remarkable spatial spillover effects [52].

LES in this study refers to the proportion of power consumption in the total energy consumption of URBs; the higher the proportion of electricity to the total household energy consumption, the lower the carbon emissions per unit building area. This work mainly analyzes the reasons from two aspects. On the one hand, power energy is a typical energy source for production and consumption in different places. It needs to be transmitted over long distances to a consumption place from the production place. The carbon emission of electric energy mainly occurs in the place of production. The use of electric energy in the place of consumption does not generate a large amount of carbon emissions. This condition shows the phenomenon of “locking of carbon emissions from electricity and energy production.” On the other hand, electricity is a typical secondary energy. Whether the power energy is clean or not is mainly determined by the primary energy used to produce electricity, which is affected by the proportion of various energy structures consumed in power generation. The electric energy that is produced totally or mainly by wind, water, and nuclear energy is clean energy. By contrast, power from coal-fired thermal power is not clean energy. According to data from the China Building Energy Conservation Association, the carbon emission coefficient of electrical energy consumption in URBs is 1.96 kg CO₂/Kgce. The carbon emission coefficient of electric energy is lower than the weighted average value of the carbon emission coefficient of URBs.

Table 7. Proportion of energy structure and carbon emission coefficient of URB in China in 2016.

	Coal	Natural Gas	Liquefied Petroleum Gas	Electricity	Heat
Energy consumption (10,000 tce)	401	3986	1720	13656	14159
Energy consumption ratio	1.18%	11.75%	5.07%	40.26%	41.74%
Carbon emission coefficient (kgCO2/kgce)	2.53	1.65	1.81	1.96	2.51
Weighted average	2.15

shows the proportion of the energy structure and carbon emission coefficient of the URBs in China in 2016.

UCL has a negative impact on the URBCE intensity, but the value is small. As stated in Section 3.1, the level of household consumption has a dual impact on the URBCE intensity. In the first stage, URBCE will increase with household consumption levels. In the second stage, after the residents’ consumption level reaches a certain level, it represents an improvement in the regional economic development level, cultural quality level, and urban construction. Residents’ energy saving awareness, life concepts, and living habits began to change; they had started pursuing an environmentally friendly lifestyle. At this time, the carbon emission intensity of residential buildings began to decrease with the improvement of residents’ living consumption levels. According to the estimated results of SDM in this work, China is in the early phase of the second stage. The elasticity coefficient of UCL on URBCE intensity is negative, but the value is small, and the negative impact is unremarkable.

6. Conclusions

In this study, the carbon emission models of residential buildings in South and North China cities are constructed, and the carbon emission intensities of these cities are calculated. The variation characteristics of URBCE intensity are analyzed by time series trends. The panel data of URBCE intensity in 30 regions of China from 2000 to 2016 are analyzed by a spatial econometric model (SDM) to discuss the factors that influence URBCE intensity further. The conclusions are presented as follows.

(1) The regions with high-energy consumption are mainly climate-affected northern district heating zones, with a large base of total energy consumption and high energy saving potential. The study found that UR and PLA are important variables affecting the energy consumption of URBs in high-energy consumption areas. It is recommended to control the total number of URBs, improve the energy efficiency of district heating in URBs, and introduce a carbon tax or carbon emission quota policy to promote the use of clean energy and reduce the total energy consumption.

(2) The consumption level of urban residents and the URB area will not only reduce the URBCE intensity in a region but also cause the URBCE intensity in adjacent areas. The level of household consumption has a dual impact on the LCI of URBs. In the first stage, URBCE will increase with household consumption levels. In the second stage, after the residents’ consumption level reaches a certain value, it represents an improvement in the regional economic development, cultural quality level, and urban construction. At this time, the carbon emission intensity of residential buildings will decrease with the improvement of residents’ living consumption levels. China is in the early phase of the second stage. UCL can weaken the carbon emission intensity, but the effect is relatively weak.

Therefore, focusing on the improvement of residents’ awareness of energy saving in areas with medium energy consumption, and the promotion and application of energy saving technologies in areas with low energy consumption may achieve better emission reduction effects.

(3) UR and PGDP have no significant impact on URBCE intensity in adjacent areas, whereas HDD changes can affect the URBCE intensity in adjacent areas. Therefore, while proposing national energy-saving and emission-reduction policies, establish measures for total energy consumption and intensity of urban residential buildings in the region. Not only can reduce the cost of policy formulation, but also achieve better energy conservation and consumption reduction effects.

This study contributes to the research on estimation methods of building carbon emissions, quantification of building climate zone indicators, and analysis of spatial differences of carbon emissions in URBs, thereby possibly providing a reference for relevant policy making. However, certain shortcomings remain. For instance, some regional data are missing in some of the covered years. Although linear interpolation is used to complete the data, its accuracy still needs to be verified further. Moreover, the spatial difference between the north heating and south non-heating areas can be more deeply investigated in follow-up studies.