Spatiotemporal Distribution of Soil Thermal Conductivity in Chinese Loess Plateau

Abstract

The Chinese Loess Plateau (CLP) is ecologically fragile, and water resources are extremely scarce. Soil thermal conductivity (λ) is a vital parameter for controlling surface heat transfer and is the key to studying the energy exchange and water balance of the soil surface. The objective of this study is to investigate the spatial distribution characteristics of soil thermal conductivity on the Loess Plateau. The research primarily employed soil heat transfer models and the Google Earth Engine (GEE) platform for remote sensing cloud computing, compares and analyzed the suitability of six models (Cambell model, Lu Yili model, Nikoosokhan model, LT model, LP1 model, and LP2 model), and utilized the selected improved model (LT model) to analyze the spatiotemporal characteristics of thermal conductivity on the CLP, examining the impacts of soil particle composition, bulk density, elevation, moisture content, and land use on thermal conductivity. The results show that the LT model is the best in the relevant evaluation indices, with a determination coefficient (R²) of 0.84, root mean square error (RMSE) of 0.18, and relative error (RE) of 0.16. Furthermore, the λ on the CLP shows an overall trend of increasing from northwest to southeast, with a lower λ between May and August and a higher one between September and October. The λ of different land use types is as follows: built-up land > cropland > grassland > forest land > barren. The bulk density (ρ_b) and altitude mainly influence λ in the CLP. The results of this study can provide a theoretical basis for studying hydrothermal variation in the CLP, model application, energy development, and land resource use.

1. Introduction

The thermal properties of soil are crucial factors in the transfer of energy and mass between the soil and the atmosphere. The thermal condition of soil also plays a significant role in soil formation and biogeochemical processes. For instance, soil temperature impacts crop root water uptake, seed germination, plant growth, and organic matter decomposition [1]. Soil thermal conductivity (λ) is the main factor that affects the soil heat transfer process and soil temperature distribution and energy change and is important in the research of thermal environment change law, which is widely used in the hydrothermal change, pollution prevention and control, and engineering construction processes, in addition to the development and utilization of shallow geothermal energy [2,3]. The soil thermal conductivity is influenced by several factors, including soil capacity, saturation (S_r), soil texture, and land use [4,5,6,7]. These different factors are coupled and interact with each other, which results in soil thermal conductivity variability in different regions.

Soil thermal conductivity is an important input parameter for numerical models [5,8], and it can be obtained through direct measurements and model predictions. Significant progress has been made in recent years in determining soil thermal conductivity, including the thermal pulse method, line source method, and hot plate method. However, the direct measurement of soil thermal conductivity requires a vast amount of time and effort [9,10]. National and international researchers have proposed several prediction models for the rapid and efficient measurement of soil thermal conductivity. Kersten [11] proposed an empirical model of soil thermal conductivity related to ρ_b, which was one of the first methods for λ prediction. However, it had limited applicability and was found to only be suitable for low water contents. Johansen [12] proposed the concept of normalized thermal conductivity (Kersten number, K_e) based on research by Kersten [11], establishing the relationship between K_e and S_r. A new thermal conductivity model was then developed on this basis. Cote and Konrad [13] proposed a new formula for the calculation of K_e that considered the effect soil texture has on λ. Lu et al. [14] further improved the Johansen [12] model by modifying the thermal conductivity of dry soil (λ_dry) and K_e calculation method. This model simulates λ in a broader range of soil water content with higher accuracy and has been strongly promoted in China [15]. Furthermore, Su et al. [16] proposed a soil thermal conductivity model that was based on basic soil physical parameters, fully considering the effects soil particle composition and organic matter content have on λ.

With practical applications, it can often be difficult to measure the input parameter (quartz content (q)) of the majority of empirical models, which further affects the prediction and study of soil thermal conductivity on a large spatial scale. As a response to this problem, Peters-Lidard et al. [17] provided recommended q values for soils with different textures. Tarnowski et al. [18] proposed an equation for estimating q from sand content (C_sand) based on the inversion of q from measured soil thermal conductivity data. Several researchers have also developed thermal conductivity models where q is ignored. Campbell [19] proposed a model that was based on ρ_b, water content (θ), and clay content (C_clay), while Lu et al. [20] developed a model that estimated thermal conductivity based on soil texture, ρ_b, and θ. Nikoosokhan [21] similarly improved the Johansen model [12] with the proposal of a new model that was based on easily measurable parameters, including C_sand, ρ_b, and S_r.

These predictive models provide invaluable methodological support for scholars who are studying the thermal conductivity regular patterns of the CLP. However, there is a scarcity of soil thermal conductivity field measurements due to the vastness of the CLP area, the large spatial heterogeneity, and difficulties in the implementation of field measurements. The available observed data is not representative of thermal conductivity distribution across the CLP, causing uncertainty in terms of the application as input parameters for numerical models. As satellite remote sensing technology has developed, the advantages of good temporal continuity and spatial homogeneity make it an essential and effective means for studying soil thermal conductivity. Google Earth Engine (GEE) is an open-source intelligent cloud platform that is dedicated to the interpretation and operation of remotely sensed images and other spatial data [22]. It provides a significant amount of satellite image data, while also using the most advanced cloud computing capabilities of Google for dealing with a variety of hotspot issues, which significantly improves efficiency and provides convenience for users [23,24,25]. It provides a vivid presentation of the topographical and geomorphological features and hydrothermal conditions of the study area while also offering the possibility of understanding the distribution of soil thermal conductivity on the CLP.

Following an evaluation and comparison of the aforementioned models, the Lu et al. model [14] is combined with the q calculation formula that was proposed by Tarnawski et al. [18] in this study for estimating the λ at different depths. The spatial distribution of λ is mapped by using the GEE cloud platform to analyze the spatial and temporal distribution characteristics of λ and its influencing factors. The aim of this study is to establish a robust approach for analyzing soil heat transfer in the CLP, with the goal of offering theoretical underpinning for regional land utilization and agricultural productivity.

2. Materials and Methods

2.1. Study Area Overview

The Loess Plateau is situated in the central part of the Yellow River Basin (33°43′~41°16′ N, 100°54′~114°33′ E) and covers an area of approximately 634,000 km². It accounts for 6.5% of the area of China and is approximately 800 km wide from north to south and 1300 km long from east to west. Elevation ranges between 95 and 5089 m above sea level, and the terrain is highest in the west and lowest in the east (Figure 1). Most higher altitude areas are located in Qinghai province in the southwest of the Loess Plateau, whereas the lower altitude areas are located primarily in the south-central part of the Loess Plateau. The administrative area covers seven provinces (autonomous regions) in Shaanxi, Shanxi, Gansu, Qinghai, Ningxia, and Henan [26,27]. The loess distribution is concentrated and deep and accounts for 70% of global loess distribution, making it the largest loess accumulation area in the world. Average rainfall on the Loess Plateau ranges from 150 to 750 mm, most of which occurs between June and September, with a highly uneven spatial and temporal distribution and mostly high-intensity rainstorms. Average annual temperature ranges from 3.6 to 14.3 °C, with a severely cold winter and warm summer. The climate is influenced by longitude and latitude, it is governed by topography, and it has the characteristics of a typical continental monsoon climate. The soil of the Loess Plateau is characterized by wind-deposited fine particles, high porosity, loose structure, and limited water retention capacity. This results in susceptibility to severe soil erosion during heavy rainfall events. Additionally, the organic matter and nutrient content of the soil are generally low [28].

2.2. Data Sources

Five types of data were selected using the GEE platform. (1) Soil particle composition data were obtained from the USDA at six standard depths (0, 10, 30, 60, 100, and 200 cm) for sand, silt, clay, and bulk density, all at a spatial resolution of 250 m. (2) Land use type data were obtained from MCD12Q1 V6, a global land cover product obtained from MODIS Terra and Aqua surface reflectance data with a spatial resolution of 500 m using supervised classification. (3) Soil water content data were obtained from the ERA5-Land dataset provided by the European Centre for Medium-Range Weather Forecasts (ECMWF), the land dataset with a spatial resolution of 11,132 m. (4) Elevation data (DEM) were obtained using DEM data from the ShuttleRadar Topography Mission (SRTM) with a resolution of 30 m (data sources are listed in

Table 1. Summary of data sources.

Dataset	Time	Resolution	Source
Soil particle composition	2018	250 m	Machine learning predictions based on global soil profile compilation by the USDA system (https://zenodo.org/records/2525663%22, accessed on 17 May 2022)
Land use	2020	500 m	United States Geological Survey (https://lpdaac.usgs.gov, accessed on 18 May 2022)
Soil water content	2000–2020	11,132 m	European Centre for Medium-Range Weather Forecasts (https://www.ecmwf.int, accessed on 19 May 2022)
Elevation	2000	30 m	SRTM database jointly provided by NASA and the Department of Defense National Mapping Agency (NIMA) (https://srtm.csi.cgiar.org, accessed on 19 May 2022)

). (5) A total of 269 measured data points on the thermal conductivity of 12 soil textures were collected from published literature (Xiong et al., 2021 [29]; Wang et al., 2012 [30]; Zeng et al., 2018 [31]; Lu et al., 2007 [14]; Tarnawski, 2009 [18]). The data of soil moisture content, soil characteristic parameters, etc., provided in these studies are input into the model for calculation to obtain the predicted values of thermal conductivity. By comparing the calculated results with the actual measured results provided in the literature, validation can be conducted.

2.3. Soil Thermal Conductivity Model

Three thermal conductivity models were selected that ignore quartz content as a parameter: the Campbell model [19], the Nikoosokhan model [32], and the Lu et al. model [20]. In addition, three methods were chosen for obtaining quartz content: (1) Peters-Lidard et al. [17] provided the recommended values for quartz content, (2) Tarnawski et al. [18] provided a formula for the estimation of quartz content, and (3) Peters-Lidard et al. [17] approximated quartz content based on sand content. The thermal conductivity models and quartz content estimation methods are combined, as shown in

Table 2. The combined models for the estimation of thermal conductivity.

Thermal Conductivity Method	Quartz Content Estimation Method	Combined Model (Abbrev.)
Campbell [19]	No
Nikoosokhan [26]	No
Lu et al. [20]	No
Lu et al. [14]	Tarnawski et al. [18]	LT model
Lu et al. [14]	Peters-Lidard et al. [17], recommend values	LP1 model
Lu et al. [14]	Peters-Lidard et al. [17], use sand content	LP2 model

.

2.3.1. Introduction to the Model

(1)

Johansen model

Johansen [12] proposed the concept of normalized thermal conductivity (K_e), which was based on the experimental results of Kersten (1949) [11]. The relationship between λ and K_e is as follows:

$$\lambda =({\lambda }_{sat}-{\lambda }_{dry})·K_e+{\lambda }_{dry}$$

(1)

where λ_sat is saturated soil thermal conductivity (W·m⁻¹·K⁻¹). K_e with S_r is associated with the following:

$$K_e=\left\{\begin{array}{c}0.7log{S}_r+1.0 0.05<S_r\le 0.1, \\ log{S}_r+1.0, S_r>0.1\end{array}\right.$$

(2)

For the λ_sat, this is a weighted average that is based on soil porosity (n), the thermal conductivity of soil constituents, and their relative content as follows:

$${\lambda }_{sat}={\lambda }_w^n{\lambda }_s^{1-n}$$

(3)

where λ_w is thermal conductivity of water, λ_w = 0.594 W·m⁻¹·K⁻¹. λ_S is the thermal conductivity of the soil solid:

$${\lambda }_{\mathrm{s}}={\lambda }_q^q{\lambda }_0^{1-q},$$

(4)

where λ_q is the thermal conductivity of the quartz content, λ_q = 7.7.0 W·m⁻¹·K⁻¹, λ₀ is the thermal conductivity of other minerals, λ₀ = 2.0 W·m⁻¹·K⁻¹ (q > 0.2), and λ₀ = 3.0 W·m⁻¹·K⁻¹ (q ≤ 0.2);

For the λ_dry, Johansen proposed the following empirical equation:

$${\lambda }_{dry}={\displaystyle \frac{0.135{\rho }_b+64.7}{2700-0.947{\rho }_b}}$$

(5)

(2)

Lu et al. (2007) [14] model

In order to improve the accuracy of the calculated λ values from the Johansen [12] model with low water content, Lu et al. [14] proposed K_e with S_r and a new expression for the exponential function with the following:

$$K_e=exp\left(\alpha \right(1-S_r^{\alpha -1.33}\left)\right)$$

(6)

where α is an empirical parameter that is determined by the sand content and takes the values of 0.96 for coarse soils with sand content greater than 40% and 0.27 for fine soils with sand content below 40%. This gives a linear relationship, predicted as λ_dry by n:

$${\lambda }_{dry}=-bn+c$$

(7)

where b and c are empirical coefficients. At 0.2 < n < 0.6, b takes the value 0.56, and c takes the value 0.51.

(3)

Lu et al. (2014) [20] model

Lu et al. [20] were aware of the effect particle composition had on λ, considered the effect of C_clay on λ, and proposed the following exponential function as a means of representing the non-linear relationship that exists between λ and θ, soil texture, and ρ_b:

$$\lambda ={\lambda }_{dry}+\mathrm{e}\mathrm{x}\mathrm{p}(\phi -{\theta }^{-\mathsf{\tau }}),\mathsf{\theta }>0$$

(8)

where τ and φ are shape factors. λ_dry uses calculations from the Lu et al. (2007) [14] model.

τ satisfies the following relationship with the C_clay particles:

$$\tau =0.67C_{clay}+0.24$$

(9)

Parameter φ has the following multiple regression equation:

$$\phi =1.97C_{sand}+1.87{\rho }_b-1.36C_{sand}{\rho }_b-0.95$$

(10)

(4)

Nikoosokhan model

Nikoosokhan [21] improved the Johansen model [12], creating a model that proposes two new sets of linear relationships between C_sand and dry density (γ_d):

$${\lambda }_{Sat}=C_{sand}+0.1{\gamma }_d$$

(11)

$${\lambda }_{dry}=0.087C_{sand}+0.019{\gamma }_d$$

(12)

The equation for the relationship between γ_d (KN·m⁻³) and ρ_b (g·m⁻³) is as follows:

γ_d = g·ρ_b

(13)

where g is the acceleration of gravity (m⁻¹·s²).

K_e is related to S_r in the equation that was proposed by the Cote and Konrad [13] model:

$$K_e={\displaystyle \frac{K{S}_r}{1+(k-1)S_r}}$$

(14)

The parameter K in Equation (14) was obtained by the following equation and can be expressed as follows:

$$K=4.4C_{sand}+0.4$$

(15)

(5)

Campbell (1985) [19] model

Campbell (1985) [19] proposed an empirical formula for the calculation:

$$\lambda =\mathrm{A}+\mathrm{B}\theta -(\mathrm{A}-\mathrm{D})\mathrm{e}\mathrm{x}\mathrm{p}[-{\left(\mathrm{C}\theta \right)}^4]$$

(16)

$$\mathrm{A}=0.65-0.78{\rho }_b+0.60{{\rho }_b}^2$$

(17)

$$\mathrm{B}=1.06{\rho }_b$$

(18)

$$\mathrm{C}=1+2.6/{C_{clay}}^{0.5}$$

(19)

$$\mathrm{D}=0.03+0.1{{\rho }_b}^2$$

(20)

(6)

Tarnawski quartz content formula

A method for the estimation of q was proposed, based on measured soil thermal conductivity and quartz content data.

$$q=0.339+0.417C_{sand}$$

(21)

(7)

Peters-Lidard (1998) [17] quartz content recommended value

Peters-Lidard [17] noted the correlation between the q and the C_sand of soil, and based on the assumption that all sand is composed of quartz, a set of relationships between q and soil texture were derived, as can be seen in

Table 3. Recommended values for quartz content in different soil textures.

Serial Number	Soil Texture	Quartz Content
1	Sand	0.92
2	Loamy sand	0.82
3	Sandy loam	0.6
4	Sandy clay loam	0.6
5	Sand clay	0.52
6	Loam	0.4
7	Clay loam	0.35
8	Silt loam	0.25
9	Clay	0.25
10	Silty clay	0.1
11	Silty clay loam	0.1
12	Silt	0.1

. The soil textures in the Loess Plateau are mainly loamy, sandy loamy, and sandy soils.

2.3.2. Statistical Analysis

As a means of quantitatively evaluating the computational accuracy of the model, the root mean square error (RMSE), normalized RMSE (NRMSE), coefficient of determination (R²), and relative error (RE) were used for assessing the conformity of the calculated and measured values of the model [33,34,35]:

$$RMSE=\sqrt{{\sum }_{i=1}^n{(O_i-S_i)}^2/n}$$

(22)

$$NRMSE=RMSE/\left(O_{max}-O_{min}\right)$$

(23)

$$R^2=1-{\sum }_i^n{(O_i-S_i)}^2/{{\sum }_i^n{(O}_i-\overline{O})}^2$$

(24)

$$RE=\sqrt{{{\sum }_{i=1}^n(O_i-S_i)}^2/{{\sum }_{i=1}^n{O}_i}^2}$$

(25)

where O_i is the actual measured value, S_i is the model calculated value, n is the sample size, and $\overline{O}$ is the sample mean.

To analyze the dominant factors affecting the thermal conductivity of Loess Plateau soil, this study used principal component analysis (PCA) to comprehensively analyze the contribution of factors such as moisture content, soil bulk density, soil texture, land use, and altitude to the variation in soil thermal conductivity. To enhance the interpretability of the principal component analysis, this study employed a varimax rotation method for further analysis.

3. Results

3.1. Characteristics of the Distribution of Physical Properties of Soils on the Loess Plateau

3.1.1. Mechanical Composition of Soil Particles

The mechanical components of the soil are the proportion of soil particles, which is the basis for studying the various physical and chemical behaviors of the soil and reflects the degree of clay-sandiness. Soil particle size affects the arrangement between the solid particles while also influencing soil structure and compactness, thereby causing thermal conductivity changes in the soil [36]. Soils that have a higher sand content and a higher proportion of quartz generally have a higher soil thermal conductivity. Therefore, soil particle composition condition is a crucial factor that affects the thermal conductivity of soil. Variations in the vertical distribution of soil particle composition on the CLP are shown in

Table 4. Statistical table of soil characteristics of Loess Plateau obtained from GEE platform.

Soil Depth (cm)	Sand Content (%)	Silt Content (%)	Clay Content (%)	Soil Bulk Density (kg·m−3)
0	44.70	39.10	16.20	1.32
10	44.21	39.01	16.78	1.34
30	43.42	38.82	17.76	1.40
60	43.22	38.58	18.20	1.44
100	44.36	38.10	17.54	1.48
200	45.63	37.23	17.14	1.48

.

The average sand content is 44.5%, which has a tendency to decrease and then increase between a depth of 0 and 200 cm, with a minimum of 43.22% at 60 cm and a maximum of 45.63% at 200 cm. Silt content averages 38.5% and decreases with depth, ranging from 39.10% to 37.23%. Clay content demonstrates a pattern of increasing and then decreasing, with the maximum value of 18.20% occurring at 60 cm.

Variations in the spatial distribution of soil particle composition on the CLP can be seen in Figure 2. Sand content (Figure 2a) is high in the north and low in the south, and the highest value is found in the northwestern part of the CLP, where sand content reaches more than 60% in northern Shaanxi and the majority of Inner Mongolia. At the same time, it is lower in both the southern and eastern parts, with less than 50%. The spatial distribution trend of silt and clay content is opposite to the spatial distribution trend of sand content and is low in the northwest and high in the southeast. Silt content (Figure 2b) is higher in the south of Gansu and Shaanxi and lower in the northwest of the CLP. However, extreme values have been observed at the northwest edge, and silt content is more than 50% in certain parts of Ningxia and Inner Mongolia. The overall clay content (Figure 2c) of the CLP is low, with an average of just 17.5%.

3.1.2. Soil Bulk Density

Soil bulk density is the mass of solids per unit volume of soil. This is an important parameter that concerns soil health and characterizes the degree of soil compaction [29]. As can be seen in Table 4, the soil bulk density of the CLP increases as soil depth increases, with an average value of 1.41 kg·m⁻³ and a range between 1.32 and 1.48 kg·m⁻³. Bulk density has a positive correlation with λ. This is because larger bulk density results in more solid-phase substances being contained in a unit volume of soil, and the soil particles have a closer arrangement. The contact area between particles increases correspondingly, ultimately leading to increased λ. As Figure 2d shows, the Loess Plateau has a greater bulk density in the northwest, but it is lower in Shanxi, Qinghai, and the Guanzhong region of Shaanxi, and the low values are distributed in strips.

3.2. Model Comparison

A total of 269 data points were collected from the literature, which covered all 12 soil texture types [18,29,30,36]. The predicted values of each λ were then obtained based on model fitting for the evaluation and validation of the thermal conductivity model for the soil. In Figure 3, the measured and predicted values of λ are compared. The predicted λ values from the Nikoosokhan model [32] can be seen to deviate significantly from the measured values. In addition, the Lu et al. model [20] was also found to be less effective. Quartz has the highest thermal conductivity of all soil mineral components at 7.7 W·m⁻¹·K⁻¹, which means that quartz content significantly affects the λ. At the same time, the Nikoosokhan model and Lu et al. models both ignore quartz content, which results in large deviations in the λ predictions and suggests quartz content is still essential when simpler empirical models are built. Hence, this study employed the Lu–Ren model to validate Tarnawski’s proposed formula, the recommended values provided by Peters-Lidard and others, and the assumption that quartz content equals three methods of estimating quartz content based on sand content. The corresponding models for these three methods are denoted as the LT model, LP1 model, and LP2 model.

The calculated values from the six models were evaluated using measured values. The resulting model evaluation index values can be seen in

Table 5. Error analysis of different models.

Model	RMSE	NRMSE	R2	RE
Campbell model	0.22	0.14	0.77	0.21
Nikoosokhan model	0.40	0.22	0.37	0.37
Lu Yili model	0.35	0.21	0.52	0.32
LT model	0.18	0.10	0.84	0.16
LP1 model	0.23	0.11	0.73	0.21
LP2 model	0.24	0.11	0.75	0.22

, the LT model showing the lowest RMSE (0.18), NRMSE (0.10), RE (0.16), and highest R² (0.84), meaning that it exhibited the best performance. The RMSE of the Campbell [37] model was found to be better than that of the Lu et al. model [20], with respective values of 0.22 and 0.35. The LP1 and LP2 models provided similar predictions with R² values of 0.73 and 0.75 and RE values of 0.21 and 0.22. This is because the recommended values that are given for quartz content are based on quartz content being equal to sand content, and the values are calculated for different soil textures.

3.3. Characterization of the Spatial Distribution of Soil Thermal Conductivity

The spatial distribution of λ on the CLP in the 0–200 cm soil layer was calculated using the LT model. As can be seen from Figure 4, the distribution characteristics of λ at different depths are generally similar, and the λ as a whole shows a trend of increasing and then decreasing from northwest to southeast, but there is also some variability. The λ shows a very small value in the northwestern edge of the CLP, and this distribution characteristic is especially significant in the depth of 0–10 cm soil layer, while no such phenomenon occurs in the 10–200 cm soil. This is because the θ in the northwestern area is small in the 0–10 cm range, and the λ is similar to that of λ_dry, while increased θ results in dry soil particles forming a water film. The contact area between the particles is increased by the appearance of the water film. This facilitates the transfer of heat, which results in the λ increasing considerably. In addition, the λ of the eastern part of the CLP is smaller, and the low-value area is similar to the spatial distribution of the ρ_b. This indicates that there is a high correlation between bulk density and thermal conductivity regarding spatial distribution. From the division of the provincial administrative regions of China, λ is higher in Qinghai and Gansu, in addition to in southern Inner Mongolia and northern Shaanxi, while it is lower in northern Inner Mongolia, southern Shaanxi, and the majority of Shanxi.

The results of descriptive statistics on λ at different soil depths can be seen in

Table 6. Thermal conductivity of soils on the Loess Plateau.

Depth of Soil Layer (cm)	Average (W·m−1·K−1)	Minimum Value (W·m−1·K−1)	Maximum Value (W·m−1·K−1)	Standard Deviation	Coefficient of Variation
0–10	1.18	0.24	1.81	0.18	15.29
10–30	1.29	0.71	1.92	0.12	9.29
30–100	1.38	0.73	2.11	0.13	9.74
100–200	1.41	0.65	2.39	0.16	11.93

, where the λ increases as soil depth increases. The average value of λ from 0–10 cm is 1.18 W·m⁻¹·K⁻¹, and the minimum value is just 0.24 W·m⁻¹·K⁻¹. The λ varies from 0.71 to 1.92 W·m⁻¹·K⁻¹ in the 10–30 cm range and has a mean value of 1.29 W·m⁻¹·K⁻¹. The range of thermal conductivity for 30–100 cm is 0.73–2.11 W·m⁻¹·K⁻¹, and this has a mean value of 1.38 W·m⁻¹·K⁻¹. The average thermal conductivity of 100–200 cm is 1.411 W·m⁻¹·K⁻¹. The coefficient of variation was used for determining the degree of variability in the λ of each layer: CV% ≤ 10% was considered weak variability, 10% < CV% < 100% was considered moderate variability, and CV% ≥ 100% was considered substantial variability. This demonstrates that λ is moderately variable for 0–10 cm and 100–200 cm but weakly variable for 10–30 cm and 30–100 cm. It should also be noted that the increase in λ from 30–200 cm is not very variable, while the increase from 0–30 cm is more variable, which means that care should be taken when measuring λ in the field to the depth where soil samples are collected.

3.4. Analysis of Temporal Variation in Soil Thermal Conductivity

Figure 5 shows the temporal variation of thermal conductivity from 2000 to 2021. When analyzing the changes in thermal conductivity of the Loess Plateau, we can divide the changes into three typical stages to better understand the underlying reasons. The first stage (2000–2002) is characterized by a period of rising thermal conductivity, increasing from 1.27 to 1.30. This growth reflects improvements in soil moisture, favorable climatic conditions, and sufficient precipitation, which contributed to increased vegetation cover, subsequently enhancing the soil’s water retention capacity and thermal conductivity [38]. The second stage (2004–2006) is a period of decline, with thermal conductivity decreasing from 1.34 W·m⁻¹·K⁻¹ in 2003 to 1.24 W·m⁻¹·K⁻¹ in 2006. This sharp decrease might be attributed to drought or adverse climatic factors post-2004, leading to drier soil and reduced moisture levels. Furthermore, improper soil usage, such as excessive farming or vegetation destruction, may have further exacerbated the decline in thermal conductivity [39]. The third stage (2007–2021) exhibits fluctuation and recovery. Thermal conductivity fluctuated and rose between 2007 and 2012, reaching 1.30 W·m⁻¹·K⁻¹ in 2012 W·m⁻¹·K⁻¹, and then remained relatively stable around 1.26 W·m⁻¹·K⁻¹. This reflects the stability and resilience of the thermal conductivity of the Loess Plateau. Changes in this stage may be associated with the improvement of environmental restoration measures, gradual recovery of precipitation, and improvement in vegetation cover. Effective land management and governance measures also contribute to mitigating soil erosion issues, playing a crucial role in controlling and optimizing thermal conductivity.

4. Discussion

4.1. Effect of Land Use Type on Thermal Conductivity

Due to its unique geographical location and the climatic characteristics it possesses, the CLP is a critical research area for soil erosion and hydrothermal changes in China. The study of spatial and temporal changes in soil thermal conductivity and land use in this region has great reference significance in terms of ecological restoration, agricultural production, and the conservation of soil and water. The spatial distribution of land use on the CLP can be seen in Figure 6. Grassland, cropland, and forest land account for a significant proportion of the land types in this area, with the size of the land area being grassland > cropland > forest land > barren > built-up land. Grassland accounts for 63.8% of the total area, which makes it the main land use type, while cropland is mainly concentrated in the Inner Mongolia Irrigation Area, southern Gansu, and the agricultural cultivation areas of southern Shaanxi and central Shanxi, and accounts for 25% of the total area. Forest land can mainly be found in central and southern Shaanxi along the Ziwu Mountains and Huanglong Mountains, in central and western Shanxi along the Luliang Mountains, and it accounts for 6.5% of the total area. Barren is concentrated mostly in northern Inner Mongolia and accounts for 2.2%, while built-up land is more scattered and accounts for 2%. Other land use types in the CLP, including shrubs, water bodies, and wetlands, account for 0.5% of the total area, which is a relatively small proportion, and it is therefore not considered in this study.

The soil thermal conductivity of different land use types (forest land, grassland, cropland, barren, and built-up land) in the CLP was studied at a depth of 0–30 cm. The temporal variation of soil thermal conductivity was analyzed (Figure 7). Soil thermal conductivity exhibits significant variations across different land use types in the Loess Plateau. The analysis reveals that built-up land and cropland demonstrate the highest average thermal conductivity (1.324 and 1.323 W·m⁻¹·K⁻¹, respectively), followed by grassland (1.303 W·m⁻¹·K⁻¹) and forest (1.270 W·m⁻¹·K⁻¹), while barren land shows the lowest value (1.196 W·m⁻¹·K⁻¹). Temporal analysis indicates distinct seasonal patterns, with peak values occurring during summer months (July-August) and minimum values in winter (December-January), showing a variation amplitude of 0.15–0.25 W·m⁻¹·K⁻¹ across all land use types. These variations are primarily attributed to land use patterns. Human activities in urban areas and agricultural lands have modified soil structure and moisture content, leading to enhanced thermal conductivity. Meanwhile, vegetation cover in grasslands and forests maintains relatively high thermal conductivity by regulating soil moisture and organic matter content. On the other hand, barren land with limited vegetation cover and moisture retention capacity exhibits the lowest thermal performance. The coefficient of variation (CV) ranges from 5.2% to 8.7%, indicating moderate temporal stability in soil thermal conductivity across different land use types. This suggests that land use patterns play a significant role in regulating soil thermal properties through their influence on soil structure, moisture conditions, and organic matter content.

4.2. Correlation Analysis of Soil Thermal Conductivity with the Main Influencing Factors

There are several factors that influence soil thermal conductivity, such as soil properties, topographic features, and land use. Pearson correlation analysis was conducted as a means of further clarifying the correlation between λ and the influencing factors. The correlation between λ and seven environmental factors (soil water content, bulk density, sand content, silt content, land use, and altitude) can be seen in

Table 7. Correlation matrix between influencing factors and soil thermal conductivity.

	Soil Water Content	Bulk Density	Sand Content	Silt Content	Land Use	Altitude	Soil Thermal Conductivity
Soil water content	1
Bulk density	−0.64 **	1
Sand content	−0.47 **	0.32 **	1
Silt content	0.42 **	−0.34 **	−0.95 **	1
Land use	−0.21 **	0.51 **	−0.45 **	0.41 **	1		.
Altitude	−0.04	−0.29 *	0.60 **	−0.46 **	−0.67 **	1
λ	0.27 **	0.53 **	0.11 *	−0.16 **	0.26 **	−0.19 **	1

Note: * indicates a significant correlation at the p < 0.05 level and ** indicates a significant correlation at the p < 0.01 level.

. λ was found to be highly significantly negatively correlated (p < 0.01) with silt content and altitude, and highly significantly positively correlated (p < 0.01) with soil water content, bulk weight, and land use. Thermal conductivity was found to be significantly positively correlated with sand content (p < 0.05). The majority of the environmental factors showed significant or highly significant correlations with each other, which suggests there is a certain amount of overlap in the information that the factors provided, so comprehensive consideration by dimensionality reduction is required in order to form new composite indicators. Principal component analysis (PCA) could satisfy this need. At the same time, PCA could clarify the contribution of each factor to the composite index, which is the main means for determining the weight of each factor.

PCA analysis (

Table 8. Results of principal component analysis of environmental factors.

Factor Loadings	Principal Components		Rotated Component Matrix
	PC1	PC2	PC1	PC2
Soil water content	0.48	−0.71	0.13	−0.84
Bulk density	−0.24	0.91	0.17	0.92
Sand content	−0.98	0.10	−0.82	0.50
Silt content	0.92	−0.14	0.78	−0.52
Land use	0.56	0.74	0.84	0.43
Altitude	−0.70	−0.51	−0.86	−0.17
Eigenvalue	2.90	2.17
Variance	48.25	36.09
Cumulative contribution rate	48.25	84.33

) found that the cumulative variance contribution of the first and second principal components reached 84.33%, which reflects most of the information on the effect the environmental factors have on soil thermal conductivity. The rotated component loadings of altitude and land use were high in the first principal component, at −0.86 and 0.84, respectively, with a variance contribution of 48.25%, which reflects the influence of topography and land use as the dominant factors.

The bulk density rotation component had the highest load of 0.92 in the second principal component. This was followed by soil water content, with a rotation component load of −0.84 and a variance contribution rate of 36.09%, mainly reflecting the influence of soil’s physical and chemical properties.

5. Conclusions

This study improved the soil thermal model in order to better simulate the spatiotemporal characteristics of soil thermal conductivity in the Loess Plateau and further analyzed the main factors influencing its spatiotemporal variation. The main conclusions are as follows:

(1)

In response to the difficulty of obtaining quartz content on a large scale in previous models for predicting soil thermal conductivity, this study compared and analyzed three methods for estimating quartz content applied to the Lu–Ren model and optimized a suitable simulation model for soil thermal conductivity in the Loess Plateau.

(2)

Overall, soil thermal conductivity shows a trend of first increasing and then decreasing from northwest to southeast and increasing with increasing soil depth. Soil thermal conductivity fluctuated greatly before 2006 and had slight fluctuations from 2007 to 2021 but, overall, tended toward stability. The ranking of soil thermal conductivity for different land use types is as follows: construction land > arable land > grassland > woodland > bare land.

(3)

Soil thermal conductivity in the Loess Plateau is significantly correlated with soil moisture content, bulk density, sand content, silt content, land use type, and altitude; among them, bulk density and altitude are major influencing factors for spatial variation in soil thermal conductivity in the Loess Plateau.