Machine Learning–Enhanced Modeling of Stress–Strain Behavior of Frozen Sandy Soil
机器学习增强的冻砂土应力-应变行为建模

Abstract  摘要

Many experiments and computational techniques have been employed to explain the mechanical properties of frozen soils. Nevertheless, due to the substantial complexity of their responses, modeling the stress–strain characteristics of frozen soils remains challenging. In this study, artificial neural networks (ANNs) were employed for modeling the mechanical behavior of frozen soil, while different testing strategies were carried out. A database covering stress–strain data from frozen sandy soil subjected to varying temperatures and confining pressures, resulting from triaxial tests, was compiled and employed to train the model. Subsequently, different artificial neural networks were trained and developed to estimate the deviatoric stress and volumetric strain, while temperature, axial strain, and confining pressure were considered as the main input variables. Based on the findings, it can be indicated that the models effectively predict the stress–strain behavior of frozen soil with a significant level of accuracy.
许多实验和计算技术已被用于解释冻土的力学特性。然而，由于其响应的高度复杂性，模拟冻土的应力-应变特性仍然具有挑战性。在本研究中，采用人工神经网络（ANNs）对冻土的力学行为进行建模，同时实施了不同的测试策略。通过三轴试验，收集了不同温度和围压下冻砂土的应力-应变数据，并建立了一个数据库用于模型训练。随后，训练和开发了不同的人工神经网络，以估算偏应力和体积应变，同时将温度、轴向应变和围压作为主要输入变量。根据研究结果，可以表明这些模型以较高的准确度有效预测了冻土的应力-应变行为。

1. Introduction  1. 引言

Building on frozen soil presents a unique challenge for cold regions, as it supports crucial infrastructure like pipelines, railways, and even buildings. However, ensuring the stability of these structures is a complex challenge [1,2]. Unlike their unfrozen counterparts, frozen soils exhibit sensitive and complicated mechanical behavior due to temperature variations. Accurately predicting how frozen soil responds to stress is paramount for safe and sustainable construction practices in cold regions. Frozen soil is a common feature across the globe, with permafrost, a permanently frozen layer, underlying a vast 21.8% of the Northern Hemisphere’s landmass [3]. This frozen ground is not uniform, though. Classifications like short-term frozen soil and seasonal frozen soil exist, highlighting the variations in freezing duration [4]. The amount of frozen soil even develops significantly during winter, with an estimated 50% of the land area experiencing frozen conditions in the coldest month [5]. Frozen soil’s inherent complexity arises from its heterogeneity, discontinuous nature, and highly non-linear response to stress and strain. These characteristics make it particularly difficult to accurately predict its mechanical behavior under real-world conditions [6,7]. This is why exploring the mechanical response of frozen soil under varying temperatures and stress conditions becomes crucial, especially for engineering projects in cold regions.
在冻土上进行建设对寒冷地区提出了独特的挑战，因为它支撑着管道、铁路甚至建筑等关键基础设施。然而，确保这些结构的稳定性是一个复杂的挑战[1, 2]。与未冻结的土壤不同，冻土由于温度变化表现出敏感且复杂的力学行为。准确预测冻土在应力下的响应对于寒冷地区安全可持续的建设实践至关重要。冻土是全球常见的特征，永久冻土层覆盖了北半球 21.8%的陆地面积[3]。然而，这种冻土并非均匀分布。存在短期冻土和季节性冻土等分类，突显了冻结持续时间的差异[4]。在冬季，冻土面积显著增加，估计在最冷的月份有 50%的陆地面积处于冻结状态[5]。冻土固有的复杂性源于其异质性、不连续性以及对应力和应变的高度非线性响应。 这些特性使得在真实条件下准确预测其力学行为变得尤为困难[6, 7]。正因如此，探究冻土在不同温度和应力条件下的力学响应变得至关重要，尤其是在寒冷地区的工程项目中。

Researchers have explored different concepts from theoretical aspects to develop constitutive models, describing the mechanical behavior of frozen soils, such as extended hypoplastic behavior [8], the rate-independent behavior of saturated frozen soils [9], the phenomenological elastoplastic damage constitutive model [10], and other research on the mechanical behavior of frozen soils, including their strength, load type, and deformation characteristics [11]. While valuable, traditional constitutive models for frozen soil behavior face several limitations that hinder their real-world application [12]. Traditional models are only effective for specific soil types, restricting their broader use in engineering practice. Model accuracy can be good for the specific data used to develop the model; however, this cannot be the case for a different type of stress [13]. The increasing mathematical complexities of these models often translate to a large number of parameters and require model calibration. The inherent complexity makes calibration of the involved parameters challenging, obstructing their practical use in engineering [12,13]. In addition, soil mechanical behavior is a variable influenced by factors like load history, time, pore water pressure, relative density, preload, pressure field, and loading rate [14].
研究人员从理论角度探索了不同的概念，以开发描述冻土力学行为的本构模型，如扩展的亚塑性行为[8]、饱和冻土的率无关行为[9]、现象学弹塑性损伤本构模型[10]，以及其他关于冻土力学行为的研究，包括其强度、荷载类型和变形特性[11]。尽管传统冻土本构模型具有重要价值，但它们在实际应用中面临诸多限制[12]。传统模型仅对特定土类有效，限制了其在工程实践中的广泛应用。模型在用于开发该模型的特定数据上可能具有良好精度；然而，对于不同类型的应力，情况并非如此[13]。这些模型日益增加的数学复杂性通常意味着需要大量参数并进行模型校准。固有的复杂性使得涉及参数的校准变得困难，阻碍了其在工程中的实际应用[12, 13]。 此外，土壤力学行为是一个受荷载历史、时间、孔隙水压力、相对密度、预压、压力场和加载速率等因素影响的变量[14]。

Machine learning-based approaches, with their ability to learn from datasets and use computing capabilities, offer an alternative to traditional mechanical behavior modeling and techniques. In the past decade, there has been a notable increase in interest in data-driven machine learning approaches within civil and geotechnical engineering [15,16,17]. Moreover, machine learning techniques are also successfully implemented to provide insight through frozen soil characteristics prediction. Artificial neural networks (ANNs) optimized with genetic algorithms were implemented for predicting unfrozen water content in frozen clay [18]. The long short-term memory (LSTM) approach and its combination with Monte Carlo dropout were employed to predict the stress–strain response of frozen soil, incorporating uncertainty quantification [7]. This growing body of research highlights a key distinction between traditional constitutive models and artificial intelligence-based (AI-based) approaches. While constitutive models are based upon mathematical equations and assumptions, these models excel at capturing the complex, non-linear relationships between stress and strain in soils through their powerful ability to analyze high-dimensional data [12,19]. AI-based methods are a promising road for predicting the mechanical behavior of frozen soil, offering a potentially more accurate and versatile approach compared to traditional methods.
基于机器学习的方法，凭借其从数据集中学习并利用计算能力的能力，为传统的力学行为建模和技术提供了替代方案。在过去的十年中，土木和岩土工程领域对数据驱动的机器学习方法的兴趣显著增加[15, 16, 17]。此外，机器学习技术也成功应用于通过冻土特性预测提供洞察。采用遗传算法优化的人工神经网络（ANNs）用于预测冻黏土中的未冻水含量[18]。长短期记忆（LSTM）方法及其与蒙特卡洛 dropout 的结合被用于预测冻土的应力-应变响应，并纳入不确定性量化[7]。这一不断增长的研究凸显了传统本构模型与基于人工智能（AI-based）方法之间的关键区别。 虽然本构模型基于数学方程和假设，但这些模型通过其分析高维数据的强大能力，在捕捉土壤中应力与应变之间复杂的非线性关系方面表现出色[12, 19]。基于 AI 的方法在预测冻土力学行为方面是一条有前景的道路，与传统方法相比，提供了潜在更准确且多功能的方法。

2. Research Significance  2. 研究意义

This research investigates the application of machine learning for enhancing the precision and applicability of artificial neural networks (ANNs) to mimic the stress–strain behavior of frozen sandy soils resulting from the experimental triaxial test. It aims to establish reliable and user-friendly models with an explainable and closed-form formulation for determining the stress–strain response of frozen sandy soil under different freezing temperatures and confining pressures. Thus, four different combinations of training and testing datasets, which are the outcomes of frozen samples under triaxial testing, and two different input feature combinations were also considered. To mitigate potential overfitting issues, a rigorous approach to model selection and evaluation was implemented. Model performance was assessed through comparisons with experimental data and unseen test data, and common error criteria were thoroughly evaluated and analyzed. Therefore, it was determined that if the developed models achieved an acceptable level of accuracy, they could be reliably used to provide a general estimation of the stress–strain behavior of frozen sandy soil, which is a key factor in general engineering design work. The aim is to overcome the lack of straightforward solutions for the mechanical behavior characterization of frozen sandy soils. The overall workflow is depicted in Figure 1.
本研究探讨了机器学习在提升人工神经网络（ANNs）精度与适用性方面的应用，旨在模拟冻砂土在实验三轴测试下的应力-应变行为。其目标是建立可靠且用户友好的模型，这些模型具备可解释的封闭式公式，用于确定冻砂土在不同冻结温度及围压下的应力-应变响应。为此，研究考虑了由三轴测试冻土样本结果得出的四种不同训练与测试数据集组合，以及两种不同的输入特征组合。为了缓解潜在的过拟合问题，采用了严格的模型选择与评估方法。通过对比实验数据与未见过测试数据，并全面评估与分析常见误差标准，对模型性能进行了深入考察。 因此，确定若所开发模型达到了可接受的精度水平，便可可靠地用于提供冻砂土应力-应变行为的一般估计，这是工程设计工作中的关键因素。目的在于克服冻砂土力学行为表征中缺乏直接解决方案的问题。整体工作流程如图 1 所示。

Figure 1. Overall research workflow for stress–strain behavior modeling of frozen sandy soil.
图 1. 冻砂土应力-应变行为建模的总体研究流程。

3. Materials and Methods  3. 材料与方法

The data extracted from the frozen soil triaxial test under −4 and −6 °C freezing temperatures and confining pressures of 0.3, 0.6, 0.8, and 1 MPa were compiled by Xu (2014) [20], as depicted in Figure 2. Standard sand was used in the triaxial compression tests. The sandy soil maximum and minimum diameters were 2.0 mm and 0.075 mm, respectively, and also, the diameter of 50% passing (D₅₀) was 0.7 mm. The total number of achieved data points was 212 points for stress and volumetric strain through eight distinct triaxial compression tests.
Xu (2014) [20] 整理了在−4 和−6°C 冻结温度及 0.3、0.6、0.8 和 1 MPa 围压下进行的冻土三轴试验数据，如图 2 所示。三轴压缩试验中使用了标准砂。砂土的最大和最小直径分别为 2.0 mm 和 0.075 mm，且 50%通过直径（D ₅₀ ）为 0.7 mm。通过八次不同的三轴压缩试验，共获得了 212 个应力与体积应变的数据点。

Figure 2. Stress–strain responses from triaxial test under −4 and −6 °C freezing temperatures, reproduced from Xu (2014) [20].
图 2. 在−4 和−6°C 冻结温度下三轴试验的应力-应变响应，引自 Xu (2014) [20]。

3.1. Data Preprocessing and Sampling Strategies
3.1. 数据预处理与采样策略

The considered input combinations for the model development involved the antecedent deviatoric stress and volumetric strains, denoted by q_(t−1) and ε_v(t−1), respectively; axial strain (ε_a); confining pressure (σ_c); and freezing temperature (T). Note that q_(t−1) and ε_v(t−1) are both a common variable of any constitutive model. However, these two can also be omitted, and another independent model can also be developed. At present, there is a wide variety of input parameters and frameworks employed in machine learning-based constitutive modeling, and there exists no established methodology or guidance for selecting these parameters and frameworks [12]. Both input combination cases were categorized as case 1 and case 2, representing the inclusion and exclusion of the antecedent target value in the models.
模型开发所考虑的输入组合包括前驱偏应力和体积应变，分别用 q _(t−1) 和ε _{v (t−1)} 表示；轴向应变（ε _a ）；围压（σ _c ）；以及冻结温度（T）。需要注意的是，q _(t−1) 和ε _{v (t−1)} 是任何本构模型的常见变量。然而，这两者也可以省略，并且可以开发另一个独立的模型。目前，基于机器学习的本构建模中采用了多种输入参数和框架，但在选择这些参数和框架方面尚无既定的方法或指导[12]。两种输入组合情况被分类为案例 1 和案例 2，分别代表模型中包含和排除前驱目标值。

Additionally, three scenarios for model testing were considered:
此外，模型测试考虑了三种情景：

• Scenario (I): Two out of eight tests results were reserved for the testing phase. This scenario was designed to simulate real-world conditions, where some test conditions might be new, and it illustrates the model’s ability to predict outcomes for the experiments that were not included in training phase;
  场景（I）：八项测试结果中有两项被保留用于测试阶段。该场景旨在模拟现实世界条件，其中某些测试条件可能是新的，并展示了模型对未包含在训练阶段的实验结果的预测能力；
• Scenario (II): The testing divisions were randomly selected segments of the strain–stress curves, each consisting of five consecutive data points. By putting aside pieces of the stress–strain curves, we aimed to test the model’s effectiveness in predicting partially unknown data, addressing conditions of incomplete data while performing laboratory experiments;
  场景（II）：测试部分为随机选取的应力-应变曲线片段，每段包含五个连续数据点。通过留出部分应力-应变曲线，我们旨在测试模型在预测部分未知数据时的有效性，以应对实验室实验过程中数据不完整的情况。
• Scenario (III): Data points were sampled from the entire dataset for the testing phase. Sampling from the full dataset allowed assessing the model’s generalization ability for a variety of data points with different pressure and temperature combinations. Two different sampling approaches were also employed in this scenario. The first approach was stratified sampling, in which an equal number of randomly selected points were chosen from each individual experiments. This process was selected to ensure that each class of data (each stress–strain curve) was proportionally represented in the testing phase and to ensure a confident testing phase for each different freezing and pressure condition. The second sampling approach was putting aside randomly selected points out of the whole data.
  场景（III）：测试阶段从整个数据集中采样数据点。从完整数据集中采样能够评估模型在不同压力和温度组合下对多种数据点的泛化能力。此场景中还采用了两种不同的采样方法。第一种方法是分层采样，即从每个单独实验中随机选择相同数量的点。选择此过程是为了确保每类数据（每条应力-应变曲线）在测试阶段按比例表示，并确保每种不同冻结和压力条件下的测试阶段具有可靠性。第二种采样方法是从整个数据中随机选择部分点作为测试集。

The testing data partition was approximately 20% for the entire developed model scenarios. Moreover, before utilizing the scenarios and selecting data for developing models, the data were scaled and normalized. In this process, the original data values were transformed to a range between 0.1 and 0.9, employing the linear relationship stated in Equation (1).
测试数据分区约占整个开发模型场景的 20%。此外，在利用这些场景并选择数据开发模型之前，数据经过了缩放和归一化处理。在此过程中，原始数据值被转换为 0.1 到 0.9 之间的范围，采用了公式(1)中所述的线性关系。

𝑋𝑆𝑐𝑎𝑙𝑒𝑑=[0.8×(𝑋−𝑋𝑚𝑖𝑛)(𝑋𝑚𝑎𝑥−𝑋𝑚𝑖𝑛)]+0.1$$X_{Scaled}=\left[0.8\times {\displaystyle \frac{(X-X_{min})}{(X_{max}-X_{min})}}\right]+0.1$$

(1)

where X is the variable, and X_min and X_max are the minimum and maximum of each variable, respectively. This ensured all features were presented on a comparable scale. Additionally, to avoid saturation of sigmoid transfer functions (used in the ANNs hidden layers), it is a common practice to scale data before using them for model development [21,22].
其中 X 为变量，X _min 和 X _max 分别为每个变量的最小值和最大值。这确保了所有特征都在可比较的尺度上呈现。此外，为避免 sigmoid 传递函数（用于 ANN 的隐藏层）的饱和，通常在模型开发之前对数据进行缩放[21, 22]。

3.2. Artificial Neural Networks (ANNs)
3.2. 人工神经网络（ANNs）

Artificial neural networks (ANNs) have emerged as powerful tools for modeling complex engineering problems [23]. Their ability to learn from existing patterns in experimental data enables them to predict future trends for unseen datasets. A typical ANN architecture consists of interconnected neurons organized into layers, including an input layer, multiple hidden layers (n), and an output layer (multilayer perceptron). The strength of each connection is represented by a weight value, and bias nodes are introduced for neurons in the hidden and output layers [24]. Multilayered feed-forward ANNs are frequently employed in civil and geotechnical engineering applications [25,26,27]. Their general architecture aligns with the basic ANN concept, featuring an input layer, n hidden layers, and an output layer. Interconnected neurons within these layers are linked in a feed-forward manner. While a sufficient network configuration is crucial for accurate predictions, practitioners must avoid overfitting by using too many neurons, as this can limit the model’s generalization capabilities [25]. Using ANN, models were developed in this study to estimate the deviatoric stress and volumetric strain that can lead to stress–strain responses. Neural Network Toolbox in MATLAB was implemented in this study for ANN modeling. Following the data division procedure recommended by Shahin et al. (2004) [25], around 75% of the data repository was dedicated to training the model. Moreover, feed-forward multilayer-based networks with a single hidden layer were utilized in this study.
人工神经网络（ANNs）已成为建模复杂工程问题的强大工具[23]。它们能够从实验数据中学习现有模式，从而预测未见数据集的未来趋势。典型的 ANN 架构由组织成层的互连神经元组成，包括输入层、多个隐藏层（n）和输出层（多层感知器）。每个连接的强度由权重值表示，并且在隐藏层和输出层的神经元中引入了偏置节点[24]。多层前馈 ANN 在土木和岩土工程应用中经常被使用[25, 26, 27]。它们的一般架构与基本 ANN 概念一致，包括输入层、n 个隐藏层和输出层。这些层中的互连神经元以前馈方式连接。虽然足够的网络配置对于准确预测至关重要，但从业者必须避免使用过多神经元导致的过拟合，因为这可能限制模型的泛化能力[25]。 本研究使用人工神经网络（ANN）开发了模型，以估计可能导致应力-应变响应的偏应力和体积应变。研究中采用了 MATLAB 中的神经网络工具箱进行 ANN 建模。根据 Shahin 等人（2004）[25]推荐的数据划分程序，约 75%的数据集被用于模型训练。此外，本研究使用了基于单隐藏层的前馈多层网络。

Artificial neural networks (ANNs) process data through interconnected layers of neurons. Each input signal is multiplied by a weight and summed with others in the same layer. This weighted sum is then fed to the next layer, where a similar process occurs. ANNs can have one or more hidden layers, forming a complex connection of weights. The most common types are feed forward and feed-forward backpropagation (FFBP). Training an ANN involves iteratively adjusting these weights to minimize the difference between the predicted and actual outputs. The activation function, often a sigmoid function, determines how the weighted sum influences the neuron’s output. Training algorithms like Levenberg–Marquardt (LM) optimize these weights by minimizing the squared error between the network’s prediction and the target value [28,29]. This iterative process allows ANNs to mimic complex relationships from data. For detailed mathematical formulations and governing equations of LM, readers are referred to the existing literature [22,30]. The schematic illustration of ANNs with a single hidden layer is presented in Figure 3.
人工神经网络（ANNs）通过相互连接的神经元层处理数据。每个输入信号会乘以一个权重，并与同一层的其他信号相加。这个加权和随后被传递到下一层，在那里进行类似的处理。ANNs 可以有一个或多个隐藏层，形成复杂的权重连接。最常见的类型是前馈和前馈反向传播（FFBP）。训练 ANN 涉及迭代调整这些权重，以最小化预测输出与实际输出之间的差异。激活函数（通常是 S 型函数）决定了加权和如何影响神经元的输出。Levenberg–Marquardt（LM）等训练算法通过最小化网络预测与目标值之间的平方误差来优化这些权重[28, 29]。这一迭代过程使 ANNs 能够从数据中模拟复杂关系。关于 LM 的详细数学公式和主导方程，读者可参考现有文献[22, 30]。图 3 展示了具有单一隐藏层的 ANNs 示意图。

Figure 3. The schematic ANNs architecture.
图 3. 人工神经网络（ANNs）架构示意图。

Hyperbolic tangent sigmoid transfer function (tansig) (y = 2/(1 + e^−2x) − 1) was considered for the hidden and linear (purelin) for the output layer. The use of one hidden layer to solve different nonlinear problems has been approved in the literature [31,32]. To evaluate the performance of the developed ANN models, mean squared error (MSE), linear correlation coefficient (R), and mean absolute percentage error (MAPE) were employed. The optimal number of neurons in the single hidden layer was determined through a trial-and-error process, evaluating models with hidden layer sizes ranging from 1 to 20 neurons. Considering the testing phase data performance, the best network configuration was selected. In the study, the deviatoric stress and volumetric strain were taken as two different targets, and efforts were made to develop distinct networks for each target.
双曲正切 S 型传递函数（tansig）（y = 2/(1 + e ^−2x ) − 1）被考虑用于隐藏层，而线性函数（purelin）用于输出层。文献[31, 32]已证实使用一个隐藏层来解决不同的非线性问题。为了评估所开发的人工神经网络（ANN）模型的性能，采用了均方误差（MSE）、线性相关系数（R）和平均绝对百分比误差（MAPE）。通过试错法确定了单隐藏层中神经元的最佳数量，评估了隐藏层大小从 1 到 20 个神经元的模型。基于测试阶段的数据表现，选择了最佳的网络配置。在本研究中，将偏应力和体积应变作为两个不同的目标，并努力为每个目标开发独立的网络。

3.3. Error Criteria  3.3. 误差准则

The performance of the optimized ANN model was assessed using well-known statistical criteria (Equations (2)–(7)) commonly employed to evaluate model performance and error. These metrics include correlation coefficient (R) (Equation (2)), coefficient of determination (R²) (Equation (3)), mean squared error (MSE) (Equation (4)), root mean squared error (RMSE) (Equation (5)), mean absolute percentage error (MAPE) (Equation (6)), and mean absolute error (MAE) (Equation (7)). To ensure a thorough comparison between models’ error values, these criteria were evaluated using the original, non-normalized target values after the back conversion of normalized data.
优化后的 ANN 模型性能通过常用的统计标准（公式(2)-(7)）进行评估，这些标准常用于衡量模型性能与误差。这些指标包括相关系数（R）（公式(2)）、决定系数（R ² ）（公式(3)）、均方误差（MSE）（公式(4)）、均方根误差（RMSE）（公式(5)）、平均绝对百分比误差（MAPE）（公式(6)）以及平均绝对误差（MAE）（公式(7)）。为确保模型间误差值的全面比较，这些标准是在归一化数据反向转换后，基于原始非归一化目标值进行评估的。

𝑅=⎡⎣⎢⎢⎢∑𝑛𝑖=1(𝐸𝑖−𝐸)(𝑃𝑖−𝑃)∑𝑛𝑖=1(𝐸𝑖−𝐸)2∑𝑛𝑖=1(𝑃𝑖−𝑃)2−−−−−−−−−−−−−−−−−−−−−−−−−√⎤⎦⎥⎥⎥$$R=\left[{\displaystyle \frac{\sum _{i=1}^n\left(E_i-\overline{E}\right)\left(P_i-\overline{P}\right)}{\sqrt{\sum _{i=1}^n{\left(E_i-\overline{E}\right)}^2\sum _{i=1}^n{\left(P_i-\overline{P}\right)}^2}}}\right]$$

(2)

𝑅2=⎡⎣⎢⎢⎢∑𝑛𝑖=1(𝐸𝑖−𝐸)(𝑃𝑖−𝑃)∑𝑛𝑖=1(𝐸𝑖−𝐸)2∑𝑛𝑖=1(𝑃𝑖−𝑃)2−−−−−−−−−−−−−−−−−−−−−−−−−√⎤⎦⎥⎥⎥2$$R^2={\left[{\displaystyle \frac{\sum _{i=1}^n\left(E_i-\overline{E}\right)\left(P_i-\overline{P}\right)}{\sqrt{\sum _{i=1}^n{\left(E_i-\overline{E}\right)}^2\sum _{i=1}^n{\left(P_i-\overline{P}\right)}^2}}}\right]}^2$$

(3)

𝑀𝑆𝐸=1𝑛∑𝑖=1𝑛(𝐸𝑖−𝑃𝑖)2$$MSE={\displaystyle \frac{1}{n}}\sum _{i=1}^n{\left(E_i-P_i\right)}^2$$

(4)

 𝑅𝑀𝑆𝐸=1𝑛∑𝑖=1𝑛(𝐸𝑖−𝑃𝑖)2−−−−−−−−−−−−−⎷$$RMSE=\sqrt{{\displaystyle \frac{1}{n}}\sum _{i=1}^n{\left(E_i-P_i\right)}^2}$$

(5)

𝑀𝐴𝑃𝐸=1𝑛∑𝑖=1𝑛|𝐸𝑖−𝑃𝑖𝐸𝑖|×100$$MAPE={\displaystyle \frac{1}{n}}\sum _{i=1}^n\left|{\displaystyle \frac{E_i-P_i}{E_i}}\right|\times 100$$

(6)

𝑀𝐴𝐸=1𝑛∑𝑖=1𝑛|𝐸𝑖−𝑃𝑖|$$MAE={\displaystyle \frac{1}{n}}\sum _{i=1}^n\left|E_i-P_i\right|$$

(7)

E_i and P_i represent the measured and estimated values for each data point, and n is the total number of data. In addition, 𝐸$\overline{E}$ and 𝑃$\overline{P}$ denote the mean measured and estimated values, respectively.
E 和 P 分别表示每个数据点的测量值和估计值，n 是数据的总数。此外， 𝐸$\overline{E}$ 和 𝑃$\overline{P}$ 分别表示测量值和估计值的平均值。

4. Results and Discussion
4. 结果与讨论

4.1. Optimized Developed Networks
4.1. 优化开发的网络

The number of neurons in the hidden layer was chosen based on the performance of the developed models, between 1 and 20 neurons. The outperforming networks are presented in Table 1. This table presents the configurations of the chosen models, along with their training and testing R, MSE, and MAPE values. It should be noted that these values were calculated using the normalized and scaled model outputs.
隐藏层中的神经元数量根据所开发模型的性能在 1 到 20 个神经元之间进行选择。表现优异的网络如表 1 所示。该表展示了所选模型的配置及其训练和测试的 R、MSE 和 MAPE 值。需要注意的是，这些值是使用归一化和缩放后的模型输出计算得出的。

Table 1. Optimized developed models performance for deviatoric stress and volumetric strain.

Model Target	Modeling Scenario		Input Case	Neuron No.	Training R-Value	Testing R-Value	Normalized Training MSE	Normalized Testing MSE	Training MAPE (%)	Testing MAPE (%)
Deviatoric Stress (q)	(I)		Case 1	3	0.99	0.99	0.0002350	0.0002455	2.98	3.14
			Case 2	4	0.99	0.99	0.0000380	0.0001053	1.24	1.54
	(II)		Case 1	12	0.99	0.99	0.0000032	0.0000246	0.27	0.55
			Case 2	12	0.99	0.99	0.0000108	0.0000617	0.58	0.92
	(III)	Stratified Sampling	Case 1	13	0.99	0.99	0.0000009	0.0000310	0.12	0.90
			Case 2	11	0.99	0.99	0.0000186	0.0000842	0.74	1.49
		Random Sampling	Case 1	12	0.99	0.99	0.0000012	0.0000351	0.16	0.73
			Case 2	14	0.99	0.99	0.0000077	0.0000378	0.43	0.95
Volumetric Strain (εv)	(I)		Case 1	2	0.99	0.99	0.0000101	0.0000061	1.17	0.99
			Case 2	3	0.99	0.99	0.0001397	0.0016639	4.05	9.02
	(II)		Case 1	10	0.99	0.99	0.0000010	0.0000031	0.34	0.42
			Case 2	9	0.99	0.99	0.0000047	0.0000073	0.85	0.86
	(III)	Stratified Sampling	Case 1	9	0.99	0.99	0.0000017	0.0000060	0.46	1.10
			Case 2	14	0.99	0.99	0.0000007	0.0000164	0.23	1.64
		Random Sampling	Case 1	14	0.99	0.99	0.0000006	0.0000106	0.27	1.42
			Case 2	14	0.99	0.99	0.0000014	0.0000037	0.42	0.76

Figure 4 and Figure 5 present the optimized model predictions (converted back to real, not scaled values) alongside the experimental data. Training points are depicted with blue markers, while testing points are shown as hollow circles. The red dashed line in Figure 4 and Figure 5 represents a perfect match between the estimated and experimental values (R = 1). This line signifies the ideal outcomes. Data points plotted closer to this line indicate better agreement between measured and predicted values. It can be seen that the points mostly lie in the vicinity of the ideal fit line, and almost all models’ correlation coefficients (R-value) are more than 0.99 (with only two decimal places shown in Table 2) for both targets of q and ε_v. Beyond that, additional criteria for model performance and error evaluation are investigated in the following sections.
图 4 和图 5 展示了优化后的模型预测值（转换回实际值，而非缩放值）与实验数据的对比。训练数据点以蓝色标记表示，而测试数据点则显示为空心圆圈。图 4 和图 5 中的红色虚线代表估计值与实验值之间的完美匹配（R = 1）。这条线象征着理想的结果。数据点越靠近这条线，表明测量值与预测值之间的一致性越好。可以看出，大多数数据点都位于理想拟合线附近，且几乎所有模型的相关系数（R 值）对于 q 和ε _v 两个目标均超过 0.99（表 2 中仅显示两位小数）。此外，模型性能及误差评估的其他标准将在后续章节中探讨。

Figure 4. Predicted deviatoric stress versus experimentally measured for different modeling scenarios.
图 4. 不同建模场景下预测的偏应力与实验测量值的对比。

Figure 5. Predicted volumetric train versus experimentally measured for different modeling scenarios.
图 5. 不同建模场景下的预测列车体积与实验测量结果对比。

Table 2. Models performances for deviatoric stress and volumetric strain on all data.

Target	Parameter	Case 1				Case 2
		Scenario (I)	Scenario (II)	Scenario (III)		Scenario (I)	Scenario (II)	Scenario (III)
				Stratified Sampling	Random Sampling			Stratified Sampling	Random Sampling
q	R	0.99	0.99	0.99	0.99	0.99	0.99	0.99	0.99
	R2	0.99	0.99	0.99	0.99	0.99	0.99	0.99	0.99
	MSE	0.00151	0.00020	0.00021	0.00024	0.00651	0.00056	0.00092	0.00040
	RMSE	0.03885	0.01417	0.01457	0.01560	0.08066	0.02371	0.03028	0.01991
	MAE	0.02781	0.00892	0.00646	0.00694	0.06442	0.01634	0.01945	0.01068
εv	R	0.99	0.99	0.99	0.99	0.99	0.99	0.99	0.99
	R2	0.99	0.99	0.99	0.99	0.98	0.99	0.99	0.99
	MSE	0.00133	0.00021	0.00039	0.00080	0.07634	0.00076	0.00062	0.00042
	RMSE	0.03646	0.01433	0.01978	0.02831	0.27629	0.02748	0.02483	0.02049
	MAE	0.02911	0.01124	0.01468	0.02118	0.18951	0.02005	0.01275	0.01265

4.2. Error Analysis  4.2. 错误分析

Figure 6 and Figure 7 illustrate the absolute errors between predicted and measured values (q and ε_v, respectively) for each data point in the training and testing sets. Additionally, Figure 8 presents the mean absolute error (MAE) for the testing data of each developed model. Figure 6 and Figure 8a reveal larger errors in predicted deviatoric stress for scenario (I) compared to other models in both cases 1 and 2. Additionally, both stratified and random sampled data point models for q estimation exhibited greater errors than scenario 2 within the same case class. Moreover, Figure 7 and Figure 8b demonstrate greater accuracy for ε_v predictions using scenario (II). Outstandingly, excluding the antecedent target value (Y_t−1) in case 2 of scenario (I) resulted in the least accurate model, which can be clearly seen in Figure 8b.
图 6 和图 7 分别展示了训练集和测试集中每个数据点的预测值与实测值（q 和ε _v ）之间的绝对误差。此外，图 8 呈现了各开发模型在测试数据上的平均绝对误差（MAE）。图 6 和图 8a 显示，在情景(I)中，相较于其他模型，案例 1 和案例 2 在预测偏应力时误差更大。同时，在相同案例类别中，无论是分层抽样还是随机抽样的数据点模型，q 的估计误差均大于情景 2。进一步地，图 7 和图 8b 表明，采用情景(II)进行ε _v 预测时，准确性更高。值得注意的是，在情景(I)的案例 2 中，若忽略先前的目标值（Y ₋₁ ），将导致模型准确性最低，这一点在图 8b 中清晰可见。

Figure 6. Comparing q models absolute errors for case 1 and 2 in different modeling scenarios.
图 6. 比较不同建模场景下案例 1 和案例 2 的 q 模型绝对误差。

Figure 7. Comparing ε_v models absolute errors for case 1 and 2 in different modeling scenarios.
图 7. 比较不同建模场景下案例 1 和案例 2 的ε _v 模型绝对误差。

Figure 8. Testing division MAE of the developed models for (a) deviatoric stress and (b) volumetric strain.
图 8. 所开发模型在(a)偏应力和(b)体积应变方面的测试划分 MAE。

4.3. Model Performance  4.3. 模型性能

The outputs of the models are depicted in Figure 9, Figure 10, Figure 11, Figure 12, Figure 13, Figure 14, Figure 15 and Figure 16 in the form of deviatoric stress and volumetric strain along with axial strain. Training and testing in these stress–strain curves are shown in filled and hollow markers, respectively.
模型的输出以偏应力和体积应变以及轴向应变的形式展示在图 9、图 10、图 11、图 12、图 13、图 14、图 15 和图 16 中。这些应力-应变曲线中的训练和测试分别以实心和空心标记表示。

Figure 9. Deviatoric stress versus axial strain results of scenario (I).
图 9. 情景(I)的偏应力与轴向应变结果。

Figure 10. Volumetric strain versus axial strain results of scenario (I).
图 10. 情景(I)的体积应变与轴向应变结果。

Figure 11. Deviatoric stress versus axial strain results of scenario (II).
图 11. 情景(II)的偏应力与轴向应变结果。

Figure 12. Volumetric strain versus axial strain results of scenario (II).
图 12. 情景(II)的体积应变与轴向应变结果。

Figure 13. Deviatoric stress versus axial strain results of scenario (III)—Stratified Sampling.
图 13. 情景（III）——分层抽样的偏应力与轴向应变结果。

Figure 14. Volumetric strain versus axial strain results of scenario (III)—Stratified Sampling.
图 14. 情景(III)——分层采样的体积应变与轴向应变结果。

Figure 15. Deviatoric stress versus axial strain results of scenario (III)—Random Sampling.
图 15. 情景（III）——随机采样的偏应力与轴向应变结果。

Figure 16. Volumetric strain versus axial strain results of scenario (III)—Random Sampling.
图 16. 场景(III)——随机抽样的体积应变与轴向应变结果。

These diagrams illustrate the estimated mechanical behavior of frozen soil under different freezing temperatures and confining pressures using ANNs. The plots for deviatoric stress–axial strain (q–ε_a) and volumetric strain–axial strain (ε_v–ε_a) generally agree with the experimental data from the training set. The ANN predictions for the testing sets are consistent with the experimental measurements, and this is particularly evident in the close agreement between the predicted and measured q–ε_a curves. However, some differences are observed in the volumetric strain–axial strain curves, especially for scenario (I). Table 2 presents the prediction performance using error and correlation evaluation indicators. In this table, the coefficient of determination R² is also used as an additional measure and complementary to the correlation coefficient, R-value. While R provides valuable insight into the strength and direction of the linear relationship between predicted and experimental values, R² is another effective criterion, as it quantifies the proportion of variance in the experimental data that the model can explain [33]. By showing that all R²-values exceed 0.99 (with only two decimal places shown in Table 2), it can be demonstrated that the ANN’s capability is not only to correlate well with the experimental values but also to capture nearly all variability in the mechanical behavior of frozen soil. This combined use of R and R² provides a more comprehensive validation of the model’s performance, as R² supports the predictive accuracy by confirming that the model captures almost all the variance in the data.
这些图表展示了使用人工神经网络（ANNs）估算的不同冻结温度和围压下冻土的力学行为。偏应力-轴向应变（q–ε _a ）和体积应变-轴向应变（ε _v –ε _a ）的曲线与训练集的实验数据基本吻合。对于测试集，ANN 的预测结果与实验测量值一致，特别是在预测和测量的 q–ε _a 曲线之间的高度一致性上表现尤为明显。然而，在体积应变-轴向应变曲线上观察到一些差异，尤其是在情景（I）中。表 2 通过误差和相关性评价指标展示了预测性能。在该表中，决定系数 R ² 也被用作相关系数 R 值的补充和额外衡量标准。虽然 R 为预测值与实验值之间线性关系的强度和方向提供了有价值的见解，但 R ² 是另一个有效的标准，因为它量化了模型能够解释的实验数据中的方差比例[33]。 通过展示所有 R ² 值均超过 0.99（表 2 中仅显示两位小数），可以证明 ANN 的能力不仅与实验值具有良好的相关性，还能捕捉冻土力学行为的几乎所有变异性。R 与 R ² 的这种组合使用提供了对模型性能更全面的验证，因为 R ² 通过确认模型几乎捕捉了数据中的所有方差来支持预测准确性。

Some of the results presented in Figure 11, Figure 12, Figure 13, Figure 14, Figure 15 and Figure 16 may seem to indicating an overfitting issue. These figures include both training and testing datasets, with predictions for the testing data specifically highlighted. A key indicator of overfitting would be a significant inconsistency between the model’s performance on training and testing data, which is not the case here. Also, the error metrics presented in Table 1 and Table 2 consistently show high accuracy for the unseen testing divisions, indicating that the model generalizes well in the introduced input span rather than overfitting. The testing strategies (partial stress–strain curve or data points) were designed to challenge the model’s ability to generalize in different conditions. The agreement between predictions on testing division and experimental data across the employed scenarios suggests that the model’s performance is not confined to the training data. Therefore, the model seems to effectively balance the accuracy and generalization without serious concern about a potential overfitting problem.
图 11、图 12、图 13、图 14、图 15 和图 16 中展示的部分结果可能看似存在过拟合问题。这些图表包含了训练和测试数据集，并特别突出了对测试数据的预测。过拟合的关键指标是模型在训练数据和测试数据上的表现存在显著不一致，但此处并非如此。此外，表 1 和表 2 中呈现的误差指标一致显示模型在未见过的测试部分具有高准确度，表明模型在引入的输入范围内具有良好的泛化能力，而非过拟合。设计的测试策略（部分应力-应变曲线或数据点）旨在挑战模型在不同条件下的泛化能力。在多种场景下，测试部分的预测与实验数据之间的一致性表明，模型的性能不仅局限于训练数据。因此，模型似乎在准确性和泛化之间实现了有效平衡，无需过度担忧潜在的过拟合问题。

In Table 2, all the data in the training and testing divisions are considered to evaluate the metrics. Additionally, all the developed models are compared in predicting q and ε_v for cases 1 and 2. Although all the models show a high level of accuracy, those with the least error criteria are highlighted in bold and shaded.
在表 2 中，考虑了训练和测试部分的所有数据以评估指标。此外，所有开发的模型在预测案例 1 和 2 的 q 和ε _v 时进行了比较。尽管所有模型都显示出较高的准确性，但误差标准最小的模型以粗体显示并加阴影。

4.4. Discussion on the Testing Strategy
4.4. 测试策略讨论

This study examined three distinct scenarios for model testing. Scenario (I) involved keeping individual experiments unseen, scenario (II) entailed setting aside sections of each stress–strain curve, and scenario (III) comprised sampling data points. Scenario (I) imposed the most rigorous condition, where the model was expected to predict an individual experiment based on its training on other physical states of experiments. In this scenario, confining pressures of 0.6 and 0.8 MPa at freezing temperatures of −4 and −6, respectively, were kept unseen as the test set. Since no previous data were presented to the model under these conditions, the lower accuracy compared to other modeling strategies can be justified. Moreover, stratified data sampling in models involving Y_t−1 led to improved accuracies compared to random point sampling.
本研究考察了模型测试的三种不同情景。情景（I）涉及保留个别实验未见，情景（II）要求预留每段应力-应变曲线的一部分，而情景（III）则包含数据点采样。情景（I）设定了最为严格的条件，模型需基于其他实验物理状态的训练来预测单个实验。在此情景下，冻结温度分别为-4℃和-6℃时的 0.6 MPa 和 0.8 MPa 围压被保留作为测试集。由于模型在这些条件下未接触过先前的数据，其准确性较其他建模策略低是可以理解的。此外，在涉及 Y ₋₁ 的模型中，分层数据采样相较于随机点采样，带来了更高的准确性提升。

4.5. Stress–Strain Formulation
4.5. 应力-应变公式

Unlike other studies that solely reported optimized models, this work provides the weights and bias values for each layer of the ANN model. This information allows readers to directly replicate the presented results using any spreadsheet program and calculate the estimated deviatoric stress and volumetric stress. The following equation expresses the relationship between the normalized input parameters (T, σ_c, ε_a and q_(t−1) or ε_v(t−1) in case 1 models) and the normalized output (q or ε_v).
与其他仅报告优化模型的研究不同，本工作提供了 ANN 模型每一层的权重和偏置值。这些信息使读者能够使用任何电子表格程序直接复现所呈现的结果，并计算估计的偏应力和体积应力。以下方程表达了归一化输入参数（T、σ _c 、ε _a 和 q _(t−1) 或ε _{v (t−1)} 在案例 1 模型中）与归一化输出（q 或ε _v ）之间的关系。

𝑌𝑛=𝑓𝑛1⎧⎩⎨𝑏0+∑𝑘=1ℎ⎡⎣⎢𝑤𝑘𝑓𝑛2⎛⎝⎜⎜𝑏ℎ𝑘+∑𝑖=1𝑚𝑤𝑖𝑘𝑋𝑖⎞⎠⎟⎟⎤⎦⎥⎫⎭⎬$$Y_n=f_{n1}\left\{b_0+\sum _{k=1}^h\left[w_k{f}_{n2}\left(b_{hk}+\sum _{i=1}^m{w}_{ik}{X}_i\right)\right]\right\}$$

(8)

Y_n represents the normalized target (q or ε_v); f_n1 and f_n2 are the linear and tansig transfer functions, respectively; h indicates neuron numbers in the hidden layer; X_i indicates the normalized values of each contributing variable; m is the number of contributing input variables; w_ik specifies the connecting weights between the ith input and kth neuron in the hidden layer; w_k is the associating weight between the kth neuron in the hidden layer and the output neuron; b_hk is the bias in the kth neuron of the hidden layer; and b₀ is the bias value in the output layer. Subsequently, w_ik, b_hk, w_k, and b₀, representing the weights and biases of the most accurate trained models with the least error in each case, were inserted into Equation (8), resulting in the formulations presented in Equations (9)–(12).
Y _n 表示归一化的目标（q 或 ε _v ）；f _{n 1} 和 f _{n 2} 分别是线性和 tansig 传递函数；h 表示隐藏层中的神经元数量；X 表示每个贡献变量的归一化值；m 是贡献输入变量的数量；w _ik 指定了第 i 个输入与隐藏层中第 k 个神经元之间的连接权重；w _k 是隐藏层中第 k 个神经元与输出神经元之间的关联权重；b _hk 是隐藏层中第 k 个神经元的偏置；b ₀ 是输出层中的偏置值。随后，将表示每种情况下误差最小的最准确训练模型的权重和偏置的 w _ik 、b _hk 、w _k 和 b ₀ 插入方程 (8)，从而得到方程 (9)–(12) 中所示的公式。

q—Scenario (III)—Stratified Sampling Model (Case 1)
q—情景（III）—分层抽样模型（案例 1）

𝑞(𝑛)=0.425𝐽1+0.295𝐽2−0.919𝐽3−0.156𝐽4+0.037𝐽5−0.875𝐽6−0.104𝐽7 −0.839𝐽8+6.017𝐽9+0.438𝐽10−1.144𝐽11−5.481𝐽12−4.409𝐽13+1.937$$\begin{gathered} q_{\left(n\right)}=0.425J_1+0.295J_2-0.919J_3-0.156J_4+0.037J_5-0.875J_6-0.104J_7 \\ \begin{array}{l}-0.839J_8+6.017J_9+0.438J_{10}-1.144J_{11}-5.481J_{12}\\ -4.409J_{13}+1.937\end{array} \end{gathered}$$

(9a)

⎡⎣⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢𝐽1𝐽2𝐽3𝐽4𝐽5𝐽6𝐽7𝐽8𝐽9𝐽10𝐽11𝐽12𝐽13⎤⎦⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥=𝑇𝑎𝑛𝑠𝑖𝑔⎛⎝⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎛⎝⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎡⎣⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢−1.209−0.6730.8791.5126.1200.256−1.767−0.1020.868−0.4250.1960.996−0.4522.1240.7530.050−1.2962.3607.897−1.626−2.122−3.751−1.502−0.062−2.715−0.0350.3253.186−0.4740.5946.3010.341−1.5220.799−1.8820.878−0.152−2.132−0.022−2.5532.0171.332−1.3815.6334.543−3.9611.0350.4041.829−6.1940.383−0.391⎤⎦⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎡⎣⎢⎢⎢⎢⎢𝑇𝜎𝑐𝜀𝑎𝑞(𝑡−1)⎤⎦⎥⎥⎥⎥⎥⎞⎠⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟+⎡⎣⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢−2.1032.892−0.877−0.863−5.350−3.341−3.281−0.7614.4281.855−5.1983.4791.510⎤⎦⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎞⎠⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟$$\left[\begin{array}{c}{J}_1\\ J_2\\ J_3\\ J_4\\ J_5\\ J_6\\ J_7\\ J_8\\ J_9\\ J_{10}\\ J_{11}\\ J_{12}\\ J_{13}\end{array}\right]=Tansig\left(\left(\left[\begin{array}{cccc}-1.209& 2.124& 0.325& -2.553\\ -0.673& 0.753& 3.186& 2.017\\ 0.879& 0.050& -0.474& 1.332\\ 1.512& -1.296& 0.594& -1.381\\ 6.120& 2.360& 6.301& 5.633\\ 0.256& 7.897& 0.341& 4.543\\ -1.767& -1.626& -1.522& -3.961\\ -0.102& -2.122& 0.799& 1.035\\ 0.868& -3.751& -1.882& 0.404\\ -0.425& -1.502& 0.878& 1.829\\ 0.196& -0.062& -0.152& -6.194\\ 0.996& -2.715& -2.132& 0.383\\ -0.452& -0.035& -0.022& -0.391\end{array}\right]\left[\begin{array}{c}T\\ {\sigma }_c\\ {\epsilon }_a\\ q_{(t-1)}\end{array}\right]\right)+\left[\begin{array}{c}-2.103\\ 2.892\\ -0.877\\ -0.863\\ -5.350\\ -3.341\\ -3.281\\ -0.761\\ 4.428\\ 1.855\\ -5.198\\ 3.479\\ 1.510\end{array}\right]\right)$$

(9b)

ε_v—Scenario (II) (Case 1)
ε _v —场景（II）（案例 1）

𝜀𝑣(𝑛)=−0.267𝐽1+0.143𝐽2−0.564𝐽3−3.619𝐽4−0.050𝐽5−0.068𝐽6−0.750𝐽7−0.103𝐽8+3.819𝐽9−0.111𝐽10−0.008$$\begin{gathered} {\epsilon }_{v\left(n\right)}=-0.267J_1+0.143J_2-0.564J_3-3.619J_4-0.050J_5-0.068J_6 \\ -0.750J_7-0.103J_8+3.819J_9-0.111J_{10}-0.008 \end{gathered}$$

(10a)

⎡⎣⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢𝐽1𝐽2𝐽3𝐽4𝐽5𝐽6𝐽7𝐽8𝐽9𝐽10⎤⎦⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥=𝑇𝑎𝑛𝑠𝑖𝑔⎛⎝⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎛⎝⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎡⎣⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢−0.766−2.772−0.2470.6284.6811.1240.022−2.6330.298−8.7162.5250.460−0.071−0.068−8.300−2.9410.0811.986−0.0560.276−0.111−0.7800.9300.217−0.4166.741−0.1197.0880.247−1.176−0.9090.262−0.2162.466−9.5650.825−0.8660.4002.2500.407⎤⎦⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎡⎣⎢⎢⎢⎢𝑇𝜎𝑐𝜀𝑎𝜀𝑣(𝑡−1)⎤⎦⎥⎥⎥⎥⎞⎠⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟+⎡⎣⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢3.0330.8921.1862.8506.1700.8190.157−3.3472.371−8.956⎤⎦⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎞⎠⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟$$\left[\begin{array}{c}{J}_1\\ J_2\\ J_3\\ J_4\\ J_5\\ J_6\\ J_7\\ J_8\\ J_9\\ J_{10}\end{array}\right]=Tansig\left(\left(\left[\begin{array}{cccc}-0.766& 2.525& -0.111& -0.909\\ -2.772& 0.460& -0.780& 0.262\\ -0.247& -0.071& 0.930& -0.216\\ 0.628& -0.068& 0.217& 2.466\\ 4.681& -8.300& -0.416& -9.565\\ 1.124& -2.941& 6.741& 0.825\\ 0.022& 0.081& -0.119& -0.866\\ -2.633& 1.986& 7.088& 0.400\\ 0.298& -0.056& 0.247& 2.250\\ -8.716& 0.276& -1.176& 0.407\end{array}\right]\left[\begin{array}{c}T\\ {\sigma }_c\\ {\epsilon }_a\\ {\epsilon }_{v(t-1)}\end{array}\right]\right)+\left[\begin{array}{c}3.033\\ 0.892\\ 1.186\\ 2.850\\ 6.170\\ 0.819\\ 0.157\\ -3.347\\ 2.371\\ -8.956\end{array}\right]\right)$$

(10b)

q—Scenario (III)—Random Sampling Model (Case 2)
q—情景（III）—随机抽样模型（案例 2）

𝑞(𝑛)=0.088𝐽1−0.192𝐽2−1.297𝐽3+0.582𝐽4−0.151𝐽5+5.004𝐽6−0.071𝐽7−0.048𝐽8−0.094𝐽9−5.072𝐽10−0.008𝐽11+0.054𝐽12+0.074𝐽13−4.767𝐽14+4.022$$\begin{gathered} q_{\left(n\right)}=0.088J_1-0.192J_2-1.297J_3+0.582J_4-0.151J_5+5.004J_6-0.071J_7 \\ \begin{array}{l}-0.048J_8-0.094J_9-5.072J_{10}-0.008J_{11}+0.054J_{12}\\ +0.074J_{13}-4.767J_{14}+4.022\end{array} \end{gathered}$$

(11a)

⎡⎣⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢𝐽1𝐽2𝐽3𝐽4𝐽5𝐽6𝐽7𝐽8𝐽9𝐽10𝐽11𝐽12𝐽13𝐽14⎤⎦⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥=𝑇𝑎𝑛𝑠𝑖𝑔⎛⎝⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎛⎝⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎡⎣⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢−4.301−0.551−0.102−1.222−14.1523.872−2.5435.6095.8090.915−39.0248.648−1.2316.8203.95912.6700.0240.562−13.460−0.495−2.968−15.141−6.656−0.372−69.6025.8264.2690.37712.0195.202−3.970−0.687−1.2672.0343.096−0.5805.5582.167111.4889.016−9.069−2.229⎤⎦⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎡⎣⎢⎢𝑇𝜎𝑐𝜀𝑎⎤⎦⎥⎥⎞⎠⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟+⎡⎣⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢12.15315.023−3.3000.9730.592−1.814−1.7721.6681.4641.234−39.4243.996−6.5756.485⎤⎦⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎞⎠⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟$$\left[\begin{array}{c}{J}_1\\ J_2\\ J_3\\ J_4\\ J_5\\ J_6\\ J_7\\ J_8\\ J_9\\ J_{10}\\ J_{11}\\ J_{12}\\ J_{13}\\ J_{14}\end{array}\right]=Tansig\left(\left(\left[\begin{array}{ccc}-4.301& 3.959& 12.019\\ -0.551& 12.670& 5.202\\ -0.102& 0.024& -3.970\\ -1.222& 0.562& -0.687\\ -14.152& -13.460& -1.267\\ 3.872& -0.495& 2.034\\ -2.543& -2.968& 3.096\\ 5.609& -15.141& -0.580\\ 5.809& -6.656& 5.558\\ 0.915& -0.372& 2.167\\ -39.024& -69.602& 111.488\\ 8.648& 5.826& 9.016\\ -1.231& 4.269& -9.069\\ 6.820& 0.377& -2.229\end{array}\right]\left[\begin{array}{c}T\\ {\sigma }_c\\ {\epsilon }_a\end{array}\right]\right)+\left[\begin{array}{c}12.153\\ 15.023\\ -3.300\\ 0.973\\ 0.592\\ -1.814\\ -1.772\\ 1.668\\ 1.464\\ 1.234\\ -39.424\\ 3.996\\ -6.575\\ 6.485\end{array}\right]\right)$$

(11b)

ε_v—Scenario (III)—Random Sampling Model (Case 2)
ε _v —场景（III）—随机抽样模型（案例 2）

𝜀𝑣(𝑛)=−0.556𝐽1+0.671𝐽2−9.911𝐽3−0.312𝐽4+1.153𝐽5−11.866𝐽6+1.333𝐽7+10.283𝐽8+11.820𝐽9−9.980𝐽10+10.106𝐽11+2.181𝐽12−0.171𝐽13−1.215𝐽14−0.368$$\begin{gathered} {\epsilon }_{v\left(n\right)}=-0.556J_1+0.671J_2-9.911J_3-0.312J_4+1.153J_5-11.866J_6+1.333J_7+10.283J_8+11.820J_9 \\ -9.980J_{10}+10.106J_{11}+2.181J_{12}-0.171J_{13}-1.215J_{14}-0.368 \end{gathered}$$

(12a)

⎡⎣⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢𝐽1𝐽2𝐽3𝐽4𝐽5𝐽6𝐽7𝐽8𝐽9𝐽10𝐽11𝐽12𝐽13𝐽14⎤⎦⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥=𝑇𝑎𝑛𝑠𝑖𝑔⎛⎝⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎛⎝⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎡⎣⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢−3.530−2.066−0.147−4.146−0.6676.062−1.5693.5643.0393.534−0.2071.198−1.5512.9368.3256.342−0.9160.171−1.726−7.343−1.086−0.904−3.1820.9060.9071.1311.1520.4013.2872.251−3.955−1.856−0.5650.3210.906−3.9060.3843.9473.927−0.171−2.646−0.903⎤⎦⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎡⎣⎢⎢𝑇𝜎𝑐𝜀𝑎⎤⎦⎥⎥⎞⎠⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟+⎡⎣⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢7.6415.582−3.2653.5891.221−2.9590.4640.481−1.729−0.123−3.851−0.021−2.8612.981⎤⎦⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎞⎠⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟$$\left[\begin{array}{c}{J}_1\\ J_2\\ J_3\\ J_4\\ J_5\\ J_6\\ J_7\\ J_8\\ J_9\\ J_{10}\\ J_{11}\\ J_{12}\\ J_{13}\\ J_{14}\end{array}\right]=Tansig\left(\left(\left[\begin{array}{ccc}-3.530& 8.325& 3.287\\ -2.066& 6.342& 2.251\\ -0.147& -0.916& -3.955\\ -4.146& 0.171& -1.856\\ -0.667& -1.726& -0.565\\ 6.062& -7.343& 0.321\\ -1.569& -1.086& 0.906\\ 3.564& -0.904& -3.906\\ 3.039& -3.182& 0.384\\ 3.534& 0.906& 3.947\\ -0.207& 0.907& 3.927\\ 1.198& 1.131& -0.171\\ -1.551& 1.152& -2.646\\ 2.936& 0.401& -0.903\end{array}\right]\left[\begin{array}{c}T\\ {\sigma }_c\\ {\epsilon }_a\end{array}\right]\right)+\left[\begin{array}{c}7.641\\ 5.582\\ -3.265\\ 3.589\\ 1.221\\ -2.959\\ 0.464\\ 0.481\\ -1.729\\ -0.123\\ -3.851\\ -0.021\\ -2.861\\ 2.981\end{array}\right]\right)$$

(12b)

5. Sensitivity Analysis  5. 敏感性分析

To assess the impact of different input variables on the predicted values of q and ε_v, a sensitivity analysis was performed on the ANN model’s results. This analysis, following Milne’s approach [34], utilized the current weights within the neural network to determine the relative influence of each input variable on the network output. Using the following equation, the percentage contribution (Q_ik) of each input variable (x_i) to the final output (q or ε_v) can be calculated by considering the weights connecting input neurons to hidden neurons (w_ij) and weights connecting hidden neurons to the output neuron (M_jk). The summation of weights connecting all N input neurons to a specific hidden neuron (j) is denoted by ∑𝑁𝑟=1|𝑤𝑟𝑗|$\sum _{r=1}^N\left|w_{rj}\right|$. The sum of Q_ik values for all input variables must always equal 100%.
为了评估不同输入变量对预测值 q 和ε _v 的影响，对 ANN 模型的结果进行了敏感性分析。该分析采用 Milne 的方法[34]，利用神经网络中的当前权重来确定每个输入变量对网络输出的相对影响。通过以下方程，可以计算每个输入变量(x)对最终输出(q 或ε _v )的百分比贡献(Q _ik )，方法是考虑连接输入神经元到隐藏神经元的权重(w _ij )以及连接隐藏神经元到输出神经元的权重(M _jk )。连接所有 N 个输入神经元到特定隐藏神经元的权重之和用 ∑𝑁𝑟=1|𝑤𝑟𝑗|$\sum _{r=1}^N\left|w_{rj}\right|$ 表示。所有输入变量的 Q _ik 值之和必须始终等于 100%。

𝑄𝑖𝑘=∑𝐿𝑗=1(𝑤𝑖𝑗∑𝑁𝑟=1|𝑤𝑟𝑗|𝑀𝑗𝑘)∑𝑁𝑖=1(∑𝐿𝑗=1|𝑤𝑖𝑗∑𝑁𝑟=1|𝑤𝑟𝑗|𝑀𝑗𝑘|)$$Q_{ik}={\displaystyle \frac{\sum _{j=1}^L\left(\frac{w_{ij}}{\sum _{r=1}^N\left|w_{rj}\right|}{M}_{jk}\right)}{\sum _{i=1}^N\left(\sum _{j=1}^L\left|\frac{w_{ij}}{\sum _{r=1}^N\left|w_{rj}\right|}{M}_{jk}\right|\right)}}$$

(13)

Figure 16 visually depicts the relative influence of each input variable on the target as determined by this analysis.
图 16 直观地展示了本分析确定的各输入变量对目标变量的相对影响。

As presented in Figure 17a,b, the Y_t−1, which is the antecedent deviatoric stress or volumetric strain, is the most influential contributing parameter, with the importance of 36.18% and 32.93%, respectively. However, in the case 2 model of deviatoric stress prediction, axial strain is the most effective input variable in the model, with 39.68%, and the least effective input is the confining pressure, with an importance of 29.2%. For the model of volumetric strain prediction in case 2, the temperature is the most influential parameter, with a 36.03% effect on the response, and the least effective parameter is the confining pressure, with a 31.56% influence on the target. It can be observed that out of all the input variables involved in the model, none of them were over- or underrated in the trained models.
如图 17a,b 所示，Y ₋₁ ，即先前的偏应力或体积应变，是最具影响力的贡献参数，其重要性分别为 36.18%和 32.93%。然而，在偏应力预测的案例 2 模型中，轴向应变是模型中最有效的输入变量，占比 39.68%，而最不有效的输入是围压，其重要性为 29.2%。对于案例 2 中体积应变预测的模型，温度是最具影响力的参数，对响应的影响为 36.03%，而最不有效的参数是围压，对目标的影响为 31.56%。可以观察到，在模型涉及的所有输入变量中，训练模型中没有对任何变量进行过高或过低的评价。

Figure 17. Relative importance of each contributing parameter in (a) Scenario (III)—Stratified Sampling Model (Case 1), (b) Scenario (II) (Case 1), (c) Scenario (III)—Random Sampling Model (Case 2) and (d) Scenario (III)—Random Sampling Model (Case 2).
图 17. 各贡献参数的相对重要性：(a) 情景(III)——分层抽样模型（案例 1），(b) 情景(II)（案例 1），(c) 情景(III)——随机抽样模型（案例 2）和(d) 情景(III)——随机抽样模型（案例 2）。

6. Conclusions  6. 结论

Traditional constitutive models for frozen soil behavior require specialized testing equipment and procedures, leading to time-consuming and costly processes. Moreover, the intricate internal structure of frozen soil, consisting of multiple phases, results in complex mechanical behaviors. To address these challenges, this study proposes a machine learning approach using artificial neural networks (ANNs) to predict the complete stress–strain response of frozen soils under various conditions. The ANN model was trained using a database from a previous experimental study on frozen soils, which included measurements at different temperatures and confining pressures.
传统的冻土行为本构模型需要专门的测试设备和程序，导致过程耗时且成本高昂。此外，冻土复杂的内部结构由多相组成，导致了复杂的力学行为。为了应对这些挑战，本研究提出了一种利用人工神经网络（ANNs）的机器学习方法，以预测冻土在各种条件下的完整应力-应变响应。该 ANN 模型使用了先前冻土实验研究中的数据库进行训练，该数据库包含了不同温度和围压下的测量数据。

Two different approaches were considered in the model development process: one involving the consideration of the antecedent (Y_t−1) target and addressing the time series problem (case 1) and the other using only independent input variables (case 2). Additionally, various testing division sampling strategies were employed to assess the impact of different sampling scenarios. Three scenarios were applied for testing data divisions: (I) keeping individual tests unseen, (II) setting aside sections of the stress–strain curves, and (III) selecting different data points from the curves. Furthermore, the effect of stratified sampling versus random sampling in the third scenario was investigated. The results showed that the choice of test division significantly affected the accuracy of the model, with scenario (I) resulting in the least accurate models. However, scenario (I) demonstrated potential in generalizing predictions for targets whose physical states were not previously introduced to the model. Most of the developed models exhibited outstanding accuracies in predicting the stress–strain behavior of frozen sandy soil, with correlation coefficients exceeding 0.99. Selected models, based on their accuracies relative to other developed models, were also presented as closed-form solutions, offering insights into the stress–strain behavior of frozen sandy soil.
在模型开发过程中考虑了两种不同的方法：一种涉及考虑前因（Y ₋₁ ）目标并解决时间序列问题（案例 1），另一种仅使用独立输入变量（案例 2）。此外，采用了多种测试划分采样策略来评估不同采样场景的影响。测试数据划分应用了三种场景：（I）保持个别测试未见，（II）留出应力-应变曲线的部分，（III）从曲线中选择不同的数据点。进一步研究了在第三种场景中分层采样与随机采样的效果。结果表明，测试划分的选择显著影响了模型的准确性，其中场景（I）产生的模型准确性最低。然而，场景（I）在泛化预测未先前引入模型物理状态的目标方面展示了潜力。 大多数已开发的模型在预测冻砂土的应力-应变行为方面表现出色，相关系数超过 0.99。根据其相对于其他开发模型的准确性，所选模型也以闭式解的形式呈现，为理解冻砂土的应力-应变行为提供了深入见解。

It was demonstrated in this study that artificial neural networks (ANNs) can efficiently model the stress–strain behavior of frozen sandy soils, which has proven to be challenging due to the complex and nonlinear nature of frozen soils under different temperature and pressure conditions. Traditional constitutive models are usually limited by extensive calibration requirements. In contrast, the ANN models developed here offer a reliable, practical approach that provides accurate predictions with a limited dataset, reducing the need for empirical models that may not fully capture frozen soil mechanics. It is worth mentioning that the current study is limited to using only one soil type, with eight experiments including 212 data points, which may seem less sufficient for resulting in robust conclusions. Therefore, more experimental tests are required to be included in future studies to further prove the robustness of the results obtained here and to provide a comprehensive intelligent model for frozen soil behavior. It should be noted that the classic constitutive models also lack adaptability across various soil types and conditions. The finding of this study highlights the potential of data-driven modeling to enhance practical engineering applications in cold regions, supporting more efficient and adaptable design solutions. Future works can address this study’s limitations by including additional experimental data from a wider range of soil compositions and testing conditions to further validate and expand the model’s applicability.
本研究表明，人工神经网络（ANNs）能够有效模拟冻砂土的应力-应变行为，而这一行为由于冻土在不同温度和压力条件下的复杂非线性特性，传统本构模型往往受限于广泛的校准需求。相比之下，本文开发的 ANN 模型提供了一种可靠、实用的方法，仅需有限数据集即可实现精确预测，减少了对可能无法全面捕捉冻土力学特性的经验模型的依赖。值得一提的是，当前研究仅限于使用一种土样，包含八项实验共计 212 个数据点，这可能在得出稳健结论方面显得不够充分。因此，未来研究需纳入更多实验测试，以进一步验证此处所得结果的稳健性，并为冻土行为提供一个全面的智能模型。值得注意的是，经典本构模型同样缺乏跨多种土类及条件的适应性。 本研究的发现凸显了数据驱动建模在提升寒冷地区实际工程应用中的潜力，支持更高效、适应性更强的设计解决方案。未来工作可通过纳入更广泛土壤成分和测试条件下的额外实验数据，解决本研究的局限性，进一步验证并扩展模型的适用性。