Energy consumption modeling and multi-parameter risk identification method for complex stratum shield construction

Abstract

The invention relates to an energy consumption modeling and multi-parameter risk identification method for complex stratum shield construction, and solves the problems that in the complex stratum shield construction, due to the high-dimensional and nonlinear coupling relation between geology and multiple parameters, the energy consumption prediction precision of an existing method is insufficient, and the energy consumption prediction accuracy is poor. The method comprises the steps that shield construction parameters and geological information are collected, geological features are extracted, and key parameters are screened. And constructing a Gaussian process regression model to predict tunneling specific energy, and optimizing a kernel function. Variable edge distribution is fitted, joint distribution is established through Copula, and an optimal model is selected. A sample is generated through Monte Carlo simulation, and a high-energy-consumption threshold value is set to recognize a high-risk working condition. And key parameter combinations are extracted to generate a probability density thermodynamic diagram, and risk early warning and operation optimization are realized. The method has the advantages that accurate prediction of shield construction energy consumption is achieved, the high-energy-consumption risk working condition is recognized from the probability level, and active early warning and operation optimization are supported.

Description

Energy consumption modeling and multi-parameter risk identification method for shield tunneling in complex strata Technical Field

[0001] This invention relates to the field of tunnel and underground engineering construction technology, and in particular to a method for energy consumption modeling and multi-parameter risk identification for shield tunneling in complex strata. Background Technology

[0002] Shield tunneling is a key technology for urban underground space development. Its excavation process is significantly affected by complex geological conditions (such as uneven geological hardness and groundwater), leading to strong nonlinearity, coupling, and time-varying characteristics in key construction parameters such as propulsion force, torque, and speed. This directly impacts construction safety, efficiency, and energy consumption. Therefore, there is an urgent need to establish a method that can accurately predict energy consumption and quantitatively identify risky operating conditions.

[0003] Existing analytical methods mainly include: 1. Empirical rule inference, which has a certain degree of adaptability but is difficult to generalize and lacks quantitative models; 2. Traditional statistical models or numerical simulation methods, the former is insufficient in characterizing nonlinearity and coupling relationships, and the latter has high computational cost and is difficult to dynamically evaluate.

[0004] The most prominent drawback is that existing methods are unable to effectively characterize the complex high-dimensional, nonlinear coupling dependency between geological information and multiple construction parameters, resulting in limited accuracy in energy consumption prediction. Furthermore, they are unable to accurately identify high-energy-consumption risk conditions and their corresponding key parameter combinations from a probabilistic perspective, thus restricting proactive early warning and parameter optimization of construction risks. Summary of the Invention

[0005] To achieve accurate prediction of energy consumption during shield tunneling and to identify high-energy-consumption risk conditions from a probabilistic perspective, supporting proactive early warning and operational optimization, this application provides an energy consumption modeling and multi-parameter risk identification method for shield tunneling in complex geological formations.

[0006] This application provides a method for energy consumption modeling and multi-parameter risk identification in shield tunneling construction in complex geological formations, employing the following technical solution:

[0007] A method for energy consumption modeling and multi-parameter risk identification in shield tunneling in complex geological formations includes:

[0008] S1: Collect construction parameters and corresponding geological description information during the shield tunneling process, and extract structured features from the geological description information to generate geological feature vectors;

[0009] S2, calculate the correlation between each construction parameter and the preset energy consumption index, and based on the correlation analysis results, select key construction parameters from the construction parameters;

[0010] S3 uses geological feature vectors and key construction parameters as inputs and tunneling specific energy as the prediction target to construct a Gaussian process regression model. The model is then optimized by combining multiple preset kernel functions to establish an energy consumption prediction model.

[0011] S4. Based on the key variables that constitute the input and output of the energy consumption prediction model, their marginal probability distributions are fitted respectively. On this basis, Copula joint distribution models with different dependency structures are constructed for comparison. Then, the optimal joint distribution model is selected according to the goodness-of-fit criterion.

[0012] S5. Monte Carlo simulation is performed using the optimal joint distribution model to generate simulated tunneling state samples. Based on the statistical distribution of tunneling specific energy in the samples, a high energy consumption threshold is set to identify high energy consumption working condition samples.

[0013] S6 extracts key construction parameter combinations from high-energy-consumption working condition samples, calculates their joint probability density and generates a probability density heatmap, identifies high-risk parameter combination intervals, and uses them for construction risk early warning and operation optimization.

[0014] By adopting the above technical solutions, the prediction accuracy of shield tunneling energy consumption (tunneling specific energy) under complex conditions has been significantly improved, achieving a leap from empirical inference to model quantification. More importantly, it defines high-energy-consumption risk conditions at a probabilistic level and can accurately locate the key construction parameter combination range that leads to the risk, thereby realizing proactive early warning of construction risks and optimization guidance for operating parameters. Attached Figure Description

[0015] Figure 1 is a flowchart illustrating an energy consumption modeling and multi-parameter risk identification method for shield tunneling in complex geological formations according to an embodiment of this application.

[0016] Figure 2 is a flowchart illustrating another embodiment of this application, which uses geological feature vectors and key construction parameters as inputs, tunneling specific energy as the prediction target, to construct a Gaussian process regression model, and optimizes the model by combining multiple preset kernel functions to establish an energy consumption prediction model. Detailed Implementation

[0017] The present application will be further described in detail below with reference to the accompanying drawings.

[0018] Referring to Figure 1, this application discloses a method for energy consumption modeling and multi-parameter risk identification in shield tunneling construction in complex geological formations, including:

[0019] S1 collects construction parameters and corresponding geological description information during the tunnel boring machine (TBM) excavation process, and extracts structured features from the geological description information to generate a geological feature vector. Construction parameters refer to various operational data collected by sensors during TBM excavation, including thrust, cutterhead torque, penetration depth, and tunneling speed, reflecting the TBM's operating status and construction efficiency. These parameters can be acquired in real time through the TBM's automated monitoring system. Geological description information refers to the descriptive text of strata in the geological survey report, including stratum type (e.g., clay layer, sand layer), lithological characteristics, etc., reflecting the geological conditions of the construction area. Geological description information is usually provided by the geological survey unit and recorded in text form.

[0020] To transform complex geological description information into structured feature vectors that can be used for modeling, this method employs the following three techniques:

[0021] Target coding: Based on historical construction data, the average tunneling parameters corresponding to each geological stratum type are statistically analyzed, and these averages are used as characteristic values ​​for that stratum. For example, statistics show that the average thrust for clay layers is 2400 kN, for silty clay layers it is 2600 kN, and for sand layers it is 3000 kN. If a shield tunnel segment contains "silty clay layer + sand layer", then the target coded thrust for that segment is (2600 + 3000) / 2 = 2800 kN. This method directly links geological description with construction parameters and is suitable for scenarios with abundant historical data.

[0022] One-Hot Encoding and PCA Dimensionality Reduction: First, the stratigraphic types are encoded using One-Hot encoding, such as clay layers as [1,0,0], silty clay layers as [0,1,0], and sand layers as [0,0,1]. For segments containing multiple stratigraphic types, these One-Hot vectors are summed and averaged to obtain a combined encoding. For example, the combined encoding for "clay layer + sand layer" is ([1,0,0] + [0,0,1]) / 2 = [0.5,0,0.5]. Subsequently, Principal Component Analysis (PCA) is used to reduce the dimensionality of the encoded high-dimensional vectors, extracting the principal components as geological feature vectors. Assuming the first two principal components extracted by PCA are [0.62, -0.13], the geological feature vector of this segment is [0.62, -0.13]. This method preserves the classification information of stratigraphic types while reducing dimensionality, thus improving the computational efficiency of the model.

[0023] Physical property statistical representation: This method directly utilizes the physical and mechanical properties of the formation (such as natural density, cohesion, and internal friction angle) for modeling. Physical property parameters for each formation are extracted from the geological database, and their mean, maximum, minimum, and coefficient of variation are calculated. For example, if a segment contains three formations with internal friction angles of 20°, 24°, and 32°, the characteristic vector of the internal friction angle for that segment is [mean 25.33°, maximum 32°, minimum 20°, coefficient of variation 0.20]. This method does not rely on tunneling parameters but is directly based on the physical properties of the formation, making it suitable for situations with relatively detailed geological information.

[0024] S2, calculate the correlation between each construction parameter and the preset energy consumption index, and based on the correlation analysis results, select key construction parameters from the construction parameters.

[0025] The correlation calculation and screening process is as follows: 1. Data preprocessing: First, the collected construction parameters and energy consumption index data are standardized to eliminate the influence of different dimensions and numerical ranges on the correlation calculation. The standardization formula is as follows: x 标准化 = Where x is the original data, μ is the mean of the data, and σ is the standard deviation of the data. 2. Correlation Calculation: Calculate the correlation coefficient between each construction parameter and the energy consumption index. Taking the Pearson correlation coefficient as an example, its calculation formula is: Where X and Y represent construction parameters and energy consumption indicators, respectively. and It is their mean. 3. Screening key parameters: Based on the absolute value of the correlation coefficient, screen the construction parameters that are highly correlated with energy consumption indicators. Usually, a threshold is set (such as 0.3), and only parameters whose absolute value of the correlation coefficient is greater than this threshold are considered key parameters.

[0026] S3. Using geological feature vectors and key construction parameters as inputs, and tunneling specific energy as the prediction target, a Gaussian process regression model is constructed. This model is then optimized by combining various preset kernel functions to establish an energy consumption prediction model. For details, please refer to steps S31 to S35.

[0027] The constructed Gaussian process regression (GPR) energy consumption prediction model is a hybrid architecture that nests physical mechanisms and data-driven approaches:

[0028] Physical Mechanism Layer: Based on the principles of tunnel boring machine mechanics, the product of thrust (TF), cutterhead torque (TQ), and penetration depth (PR) is used as the physical baseline for energy consumption. The formula is as follows:

[0029] Where k1 is the formation hardness correction coefficient and k2 is the mechanical efficiency coefficient, both of which are calibrated using historical data.

[0030] Data-driven layer: Employs an additive combination of the squared exponential (SE) kernel function and the Matérn3 / 2 kernel function, with the following structure: .in, , , , The hyperparameters to be trained control the global smoothness and the ability to capture local fluctuations, respectively.

[0031] Model topology: The input layer has a dimension of d (d = geological feature dimension + key construction parameter dimension), the hidden layer is induced to be mapped to an infinite-dimensional feature space through a kernel function, and the output layer is the predicted tunneling specific energy value and its 95% confidence interval.

[0032] The training data came from 5,000 rings of historical tunneling data collected by the tunnel boring machine's sensors, with a sampling frequency of 1Hz. The training process employed alternating iterative methods of maximum likelihood estimation (MLE) and Bayesian optimization.

[0033] 1. Parameter initialization: σ 2 =1.0, l=1.0, =0.5,l m =0.8, noise variance =0.1.

[0034] 2. Optimization objective: Maximize the log-likelihood function. The L-BFGS-B algorithm is used to iterate to the gradient norm. .

[0035] 3. Hyperparameter tuning: Bayesian optimization search σ 2 ∈[0.1,10], l∈[0.5,5], with the goal of minimizing the RMSE of the validation set, converges after 30 iterations.

[0036] S4. Based on the key variables that constitute the input and output of the energy consumption prediction model, their marginal probability distributions are fitted respectively. On this basis, Copula joint distribution models with different dependency structures are constructed for comparison. Then, the optimal joint distribution model is selected according to the goodness-of-fit criterion.

[0037] The specific process of fitting the marginal probability distributions of the key variables constituting the input and output of the energy consumption prediction model can be found in steps S4a to S4d. For comparing the construction of Copula joint distribution models with different dependency structures based on this, and then selecting the optimal joint distribution model according to the goodness-of-fit criterion, please refer to steps S4e to S4h, which will not be elaborated here.

[0038] S5 uses the optimal joint distribution model to perform Monte Carlo simulation, generates simulated tunneling state samples, and sets a high energy consumption threshold based on the statistical distribution of tunneling specific energy in the samples to identify high energy consumption working condition samples.

[0039] The Monte Carlo simulation and high-energy consumption identification process is as follows: 1. Generate simulation samples: Using an optimal joint distribution model (such as the VineCopula model), a large number of random samples are generated using the Monte Carlo method. Each sample contains key parameters of shield tunneling (such as thrust, cutterhead torque, penetration depth, etc.) and energy consumption indicators (such as tunneling specific energy, SEC). It is recommended that the total number of generated samples be no less than 10,000 to ensure the stability of the sample distribution and the reliability of the analysis. 2. Set a high-energy consumption threshold: Perform statistical analysis on the tunneling specific energy (SEC) in the generated simulation samples and calculate its 95th quantile (or select other high quantiles according to engineering requirements). Use this quantile value as the high-energy consumption threshold. For example, if the 95th quantile of SEC in the simulation samples is 1.2 kWh / m³, then set 1.2 kWh / m³ as the high-energy consumption threshold. 3. Identify high-energy consumption case samples: Select samples from all simulation samples whose tunneling specific energy exceeds the high-energy consumption threshold; these samples are the high-energy consumption case samples. These samples represent unfavorable conditions with high energy consumption during tunnel boring machine (TBM) construction and can be used for subsequent risk analysis and construction optimization.

[0040] S6 extracts key construction parameter combinations from high-energy-consumption operating condition samples, calculates their joint probability density and generates a probability density heatmap to identify high-risk parameter combination intervals for construction risk early warning and operation optimization. The specific process can be found in steps S6a to S6e, and will not be elaborated here.

[0041] Referring to Figure 2, a Gaussian process regression model is constructed using geological feature vectors and key construction parameters as inputs and tunneling specific energy as the prediction target. This model is then optimized by combining various preset kernel functions to establish an energy consumption prediction model, including:

[0042] S31. Select a variety of basic kernel functions from the preset kernel function set. The basic kernel functions include at least the radial basis function kernel, the Matern kernel, the rational quadratic kernel, and the white noise kernel.

[0043] The necessary steps are as follows: 1. Determine the pre-set kernel function set: In Gaussian process regression (GPR) models, the choice of kernel function is crucial to the model's performance. The kernel function defines the similarity between input data points, thus affecting the model's ability to fit and generalize the data. To construct a GPR model capable of predicting energy consumption during shield tunneling construction in complex geological formations, it is necessary to select several basic kernel functions from the pre-set kernel function set. These basic kernel functions should include at least the following:

[0044] Radial Basis Function Kernel (RBF Kernel): Used to capture smooth variations in data, suitable for continuous and smooth functional relationships.

[0045] Matérn kernel: Used to capture non-smooth changes in data, offering greater flexibility and suitable for data with discontinuous or non-smooth relationships.

[0046] Rational quadratic kernel (RQ kernel): It combines the characteristics of RBF kernel and polynomial kernel, and can capture more complex nonlinear relationships.

[0047] White noise kernel: Used to capture the noisy part of the data. It is usually used in combination with other kernel functions to improve the robustness of the model to noise.

[0048] 2. Kernel function parameter initialization:

[0049] Each kernel function has its specific parameters, which need to be optimized during model training. In step S31, these parameters need to be initialized to ensure smooth subsequent model training. The specific parameter initialization process is as follows:

[0050] Radial basis function kernel (RBF kernel):

[0051] Parameters: length scale l and signal variance .

[0052] Initialization: Typically, the length scale l is initialized to 1.0, and the signal variance... Initialized to 1.0.

[0053] Matérn core:

[0054] Parameters: smoothing parameter ν, length scale l, and signal variance .

[0055] Initialization: The smoothing parameter ν can be chosen from common values ​​(such as 1.5 or 2.5), the length scale l is initialized to 1.0, and the signal variance... Initialized to 1.0.

[0056] Rational quadratic kernel (RQ kernel):

[0057] Parameters: shape parameter α, length scale l, and signal variance .

[0058] Initialization: Shape parameter α is initialized to 1.0, length scale l is initialized to 1.0, and signal variance is... Initialized to 1.0.

[0059] White noise kernel:

[0060] Parameter: Noise variance .

[0061] Initialization: Noise Variance Initialize to 0.1 (a smaller value, indicating a lower noise level).

[0062] S32, combine the basic kernel functions in pairs to generate multiple hybrid kernel function configurations. The specific process can be found in steps S321 to S322, and will not be elaborated here.

[0063] S33 uses geological feature vectors and key construction parameters as inputs, and tunneling specific energy as the prediction target. It constructs and trains corresponding Gaussian process regression sub-models using various hybrid kernel function configurations.

[0064] Gaussian process regression (GPR) is a probabilistic regression method based on Gaussian processes. It uses kernel functions to capture the nonlinear relationship between input features and output variables, providing predicted values ​​and estimates of their uncertainties. Hybrid kernel function configurations are combinations of basic kernel functions (such as RBF, Matérn, rational quadratic kernels, etc.) through addition or multiplication, used to enhance the model's fitting ability.

[0065] Input features:

[0066] Geological feature vector: Structured geological features generated through step S1.

[0067] Key construction parameters: Construction parameters that are highly related to energy consumption, selected through step S2.

[0068] Output target: Tunneling specific energy (SEC), as the target variable for energy consumption prediction.

[0069] The necessary processes are described below: 1. Data Preparation: Combine geological feature vectors and key construction parameters into an input feature matrix X, and use tunneling specific energy (SEC) as the target variable y. 2. Model Initialization and Training: For each hybrid kernel function configuration, initialize the kernel function parameters, construct the GPR model, and optimize the parameters using maximum likelihood estimation (MLE). For example, for the RBF + white noise kernel:

[0070] Kernel function formula: .

[0071] Parameter initialization: , , .

[0072] Training process: Optimize the parameters by maximizing the log-likelihood function to obtain the optimized parameters (e.g., ...). =0.9, =1.2, =0.08).

[0073] 3. Repeat other configurations: Repeat the above process for other hybrid kernel configurations (such as Matérn + rational quadratic kernel, RBF + Matérn kernel) to optimize kernel function parameters.

[0074] 4. Model storage: Saves the optimized model and its parameters for each configuration to prepare for subsequent performance evaluation.

[0075] S34. Evaluate the predictive performance of each sub-model based on the validation dataset. The evaluation metrics for predictive performance include the coefficient of determination and the root mean square error.

[0076] The necessary process is as follows: 2.1 Data preparation: Divide the verification dataset from the shield tunneling construction data, which includes input features (geological feature vectors and key construction parameters) and target variables (tunneling energy SEC).

[0077] 2.2 Model Prediction:

[0078] For each GPR sub-model under the hybrid kernel function configuration:

[0079] 1. Use the trained model to predict the input features of the validation dataset to obtain the predicted tunneling ratio. .

[0080] 2. Record the actual tunneling energy y.

[0081] 2.3 Calculation of evaluation indicators:

[0082] Calculate the coefficient of determination (R²) for each sub-model. 2 ) and root mean square error (RMSE).

[0083] Coefficient of determination (R) 2 ): Among them, y i This is the actual value. It is a predicted value. It is the average of the actual values.

[0084] Root Mean Square Error (RMSE): .

[0085] 2.4 Performance Evaluation and Recording:

[0086] According to R 2 The RMSE values ​​are used to evaluate the predictive performance of each sub-model, and the results are recorded. Typically, R is chosen. 2 Models with higher RMSE and lower RMSE are considered to have better performance.

[0087] S35. Based on the evaluation index of prediction performance, the Gaussian process regression sub-model with the best prediction performance is selected from all sub-models as the final energy consumption prediction model.

[0088] The necessary procedures are as follows:

[0089] 2.1 Collect performance metrics for each sub-model: Obtain the R-value of each GPR sub-model from step S34. 2 The RMSE values ​​were compiled into a table for comparison.

[0090] 2.2 Selecting the optimal model:

[0091] The optimal model is selected according to the following rules:

[0092] Prefer R 2 The model with the highest value.

[0093] In R 2 In similar cases, choose the model with the smallest RMSE value.

[0094] 2.3 Confirm and record the optimal model: After determining the optimal model, record its kernel function configuration, parameters and performance indicators as the final energy consumption prediction model.

[0095] Combining basic kernel functions in pairs to generate multiple hybrid kernel function configurations includes:

[0096] S321 generates all possible unordered pairwise combinations based on the selected basic kernel functions;

[0097] The process is as follows:

[0098] 1. Determine the set of basic kernel functions: Assume that the RBF kernel, Matérn kernel, rational quadratic kernel (RQ), and white noise kernel are selected.

[0099] 2. Generate Combinatorial Pairs: From these kernel functions, generate all possible pairwise combinations, regardless of order. Specific combinations are as follows:

[0100] RBF+Matérn;

[0101] RBF+RQ;

[0102] RBF+white noise;

[0103] Matérn+RQ;

[0104] Matérn+white noise;

[0105] RQ+white noise;

[0106] Result: The generated unordered pairwise combinations will be used in the next step of generating hybrid kernel function configurations, enhancing the flexibility and adaptability of the model.

[0107] S322 generates a corresponding hybrid kernel configuration for each pair of unordered combinations by performing kernel function addition or kernel function multiplication.

[0108] The process is as follows:

[0109] Kernel function addition: This operation adds two kernel functions together to capture two different characteristics simultaneously. For example, for the RBF and Matérn kernels:

[0110] .

[0111] Kernel multiplication: Multiplying two kernel functions enhances nonlinear characteristics. For example, for RBF and rational quadratic kernels (RQ): .

[0112] Based on the key variables constituting the input and output of the energy consumption prediction model, their marginal probability distributions are fitted respectively, including:

[0113] S4a, which determines the set of key variables for the marginal probability distribution to be fitted. This set includes the geological feature vector and key construction parameters that constitute the input of the energy consumption prediction model, as well as the tunneling specific energy that constitutes its output.

[0114] The process is as follows:

[0115] 1. Input variables:

[0116] Geological feature vector: Structured geological features generated through step S1, including target encoding, PCA dimensionality reduction results, and physical property statistical features.

[0117] Key construction parameters: Construction parameters that are highly related to energy consumption and selected through step S2, such as thrust (TF), cutterhead torque (TQ), and penetration (PR).

[0118] 2. Output variables:

[0119] Specific Energy Per Mining (SEC): The energy consumed by a tunnel boring machine per unit volume of soil, used as the output target of the energy consumption prediction model.

[0120] S4b represents each variable in the set of key variables. Multiple candidate distributions are selected from a predefined family of candidate probability distributions for fitting. The family of candidate probability distributions includes various distributions such as normal distribution, log-normal distribution, Weibull distribution, Beta distribution, Logistic distribution, and Gumbel distribution.

[0121] The process is as follows:

[0122] 1. Determine the set of key variables: Based on step S4a, the set of key variables includes each component of the geological feature vector, key construction parameters, and tunneling specific energy (SEC).

[0123] 2. Preset candidate probability distribution families, and select the following common probability distributions as candidates:

[0124] Normal distribution;

[0125] Log-Normal distribution;

[0126] Weibull distribution;

[0127] Beta distribution;

[0128] Logistic distribution;

[0129] Gumbel distribution;

[0130] 3. Fit candidate distributions for each variable:

[0131] For each variable in the set, fit the candidate distribution described above.

[0132] The fitting process typically involves estimating the parameters of the distribution so that the distribution best matches the actual data of the variable.

[0133] S4c, based on at least one of the Akaike information criterion, Bayesian information criterion and log-likelihood value, evaluate the fitting effect of each candidate distribution fitted for each variable in step S4b.

[0134] The necessary procedures are as follows:

[0135] 1. Calculate the log-likelihood value: The log-likelihood value measures the probability density of data under a given distribution parameter. The higher the value, the better the fit.

[0136] For each candidate distribution of each variable, calculate its log-likelihood value:

[0137] .in, It is the probability density function under parameter θ, x i These are data points.

[0138] 2. Calculate the Akaike Information Criterion (AIC):

[0139] AIC is used for model selection, taking into account the complexity of the model (number of parameters). The smaller the value, the better the model.

[0140] Calculation formula: Where k is the number of model parameters.

[0141] 3. Calculate the Bayesian Information Criterion (BIC):

[0142] BIC is also used for model selection, and it penalizes model complexity more strictly than AIC. The smaller the value, the better the model.

[0143] Calculation formula: Where n is the number of data points.

[0144] 4. Evaluate the fitting effect:

[0145] For each candidate distribution of each variable, calculate the above goodness-of-fit index.

[0146] Comparing these indicators, the distribution with the highest log-likelihood value and the smallest AIC and BIC is selected as the optimal distribution.

[0147] S4d selects the optimal marginal probability distribution from multiple candidate distributions fitted to each variable based on the goodness-of-fit index.

[0148] The process is as follows:

[0149] 1. Review of goodness-of-fit indices:

[0150] Log-likelihood value: The higher the value, the better the fit.

[0151] Akaike Information Criterion (AIC): The smaller the value, the better the model.

[0152] Bayesian Information Criterion (BIC): The smaller the value, the better the model.

[0153] 2. Select the optimal distribution:

[0154] For each variable, compare the goodness-of-fit indices of all its candidate distributions.

[0155] The distribution with the highest log-likelihood value is usually preferred.

[0156] Given similar log-likelihood values, choose the distribution with the smallest AIC and BIC.

[0157] For each variable in the set of key variables, multiple candidate distributions are selected from a predefined family of candidate probability distributions for fitting, including:

[0158] S4b1, based on the data distribution characteristics of each variable in the set of key variables, dynamically match or select one or more corresponding candidate probability distribution types for each variable from the preset candidate probability distribution family; the data distribution characteristics are determined at least by skewness, kurtosis and empirical distribution shape.

[0159] The necessary procedures are as follows:

[0160] 1. Determine the set of key variables:

[0161] This includes the components of the geological feature vector, key construction parameters, and tunneling specific energy (SEC).

[0162] 2. Analyze the data distribution characteristics:

[0163] Skewness: Measures the asymmetry of data distribution. Positive values ​​indicate right skewness, and negative values ​​indicate left skewness.

[0164] Kurtosis: Measures the steepness or flatness of a data distribution. Kurtosis indicates that the data has a sharper peak and a thicker tail.

[0165] Empirical distribution pattern: Observe the overall pattern of the data by drawing histograms or empirical cumulative distribution functions (ECDF).

[0166] 3. Dynamic matching of candidate distribution:

[0167] Based on skewness, kurtosis, and empirical distribution shape, select one or more suitable candidate probability distribution types for each variable.

[0168] Normal distribution: suitable for symmetrical and unbiased data.

[0169] Log-normal distribution: Applicable to right-skewed data.

[0170] Weibull distribution: suitable for data with long or short tails.

[0171] Beta distribution: suitable for data that varies within the range [0,1].

[0172] Logistic distribution: suitable for symmetrical data but with thicker tails than the normal distribution.

[0173] Gumbel distribution: suitable for extreme value data.

[0174] S4b2, for each candidate probability distribution type selected for each variable, uses the maximum likelihood estimation method or the method of moments estimation to estimate parameters based on the historical observation data of that variable, and completes the fitting of the candidate distribution.

[0175] The process is as follows:

[0176] 1. Determine the candidate distribution and its parameter form:

[0177] According to step S4b1, multiple candidate distribution types (such as normal distribution, log-normal distribution, etc.) have been selected for each variable.

[0178] Each distribution has its specific parameter form, for example:

[0179] Normal distribution: parameters are mean (μ) and standard deviation (σ).

[0180] Log-normal distribution: The parameters are the logarithmic mean (μ) and the logarithmic standard deviation (σ).

[0181] Weibull distribution: The parameters are shape parameter (k) and scale parameter (λ).

[0182] 2. Selecting a parameter estimation method:

[0183] Maximum Likelihood Estimation (MLE): This method estimates distribution parameters by maximizing the likelihood function and is applicable to most distribution types.

[0184] Method of Moments (MoM): This method estimates parameters by matching sample moments (such as mean and variance) with the distribution moments. It is applicable to certain specific distributions.

[0185] 3. Parameter estimation and fitting:

[0186] For each candidate distribution of each variable, the distribution parameters are calculated using the selected parameter estimation method.

[0187] Use the calculated parameters to fit the candidate distribution.

[0188] Based on this, Copula joint distribution models with different dependency structures were constructed and compared. Then, the optimal joint distribution model was selected according to the goodness-of-fit criterion, including:

[0189] S4e constructs candidate Copula joint distribution models with at least two different dependency structures in parallel, wherein the dependency structure type is selected from a set consisting of C-vine structure, D-vine structure and R-vine structure.

[0190] The specific process is as follows:

[0191] 1. Determine the dependency structure type:

[0192] C-vine structure: Centered on a core variable, the dependencies of other variables revolve around this core variable, suitable for scenarios with a clear central variable.

[0193] D-vine structure: There are sequential dependencies between variables, which is suitable for scenarios where the dependency relationship between variables has a clear order.

[0194] R-vine structure: more flexible, allows correlations between arbitrary variables, and is suitable for complex dependencies.

[0195] 2. Construct candidate models:

[0196] Based on the number of key variables and the characteristics of dependencies, select at least two different dependency structure types.

[0197] For each dependency structure, construct the corresponding Copula joint distribution model.

[0198] S4f: For each candidate Copula joint distribution model constructed in step S4e, based on the optimal marginal probability distribution selected for each key variable in step S4d, the maximum likelihood estimation method is used to solve for the Pair-Copula function parameters of the model to complete the model parameter estimation.

[0199] The process is as follows:

[0200] 1. Review of candidate models: Based on step S4e, candidate Copula joint distribution models with different dependency structures (such as C-vine, D-vine, R-vine) have been constructed.

[0201] 2. Review the optimal marginal probability distribution: According to step S4d, the optimal marginal probability distribution for each key variable has been determined.

[0202] 3. Construct a joint distribution model:

[0203] For each candidate Copula model, the optimal marginal probability distribution of the key variables is transformed into a uniform distribution (using the cumulative distribution function CDF).

[0204] The joint distribution model is constructed step by step using the Pair-Copula function. The Pair-Copula function describes the dependency between two variables.

[0205] 4. Maximum Likelihood Estimation (MLE):

[0206] Construct a likelihood function and estimate the parameters of the Pair-Copula function by maximizing the likelihood function.

[0207] For each Pair-Copula function, maximize the following likelihood function:

[0208] .

[0209] in, It is a Pair-Copula function, u i and v i The transformed uniformly distributed data is θ, which is the parameter to be estimated.

[0210] 5. Parameter estimation:

[0211] Numerical optimization methods (such as Newton's method or quasi-Newton's method) are used to solve for the parameter θ of the likelihood function.

[0212] Repeat the above process for each Pair-Copula function to complete the parameter estimation of the entire Copula joint distribution model.

[0213] S4g, based on a unified set of goodness-of-fit criteria, quantifies and evaluates the overall goodness of fit of each candidate Copula joint distribution model for which parameter estimation has been completed in step S4f. The set of goodness-of-fit criteria includes the Akaike information criterion, the Bayesian information criterion, and the log-likelihood value.

[0214] The process is as follows:

[0215] 1. Review the models whose parameters have been estimated: In step S4f, the Pair-Copula function parameters of each candidate Copula joint distribution model have been estimated.

[0216] 2. Calculate the log-likelihood value: The log-likelihood value measures how well the model fits the data; a higher value indicates a better fit. For each candidate model, calculate the overall log-likelihood value:

[0217] .in, It is the joint probability density function under parameter θ, x i It is observational data.

[0218] 3. Calculate the Akaike Information Criterion (AIC): AIC is used for model selection. It takes into account the complexity of the model (number of parameters). The smaller the value, the better the model.

[0219] Calculation formula: Where k is the total number of model parameters.

[0220] 4. Calculate the Bayesian Information Criterion (BIC): BIC is also used for model selection. It penalizes model complexity more strictly than AIC. The smaller the value, the better the model.

[0221] Calculation formula: Where n is the number of data points.

[0222] 5. Evaluate the overall goodness of fit:

[0223] For each candidate model, calculate the goodness-of-fit index mentioned above.

[0224] Record the log-likelihood value, AIC, and BIC for each model.

[0225] S4h: Based on the preset quantitative comparison rules, compare the evaluation results of each candidate Copula joint distribution model according to the goodness-of-fit criterion set obtained in step S4g, select the candidate model with the best comprehensive evaluation result, and determine it as the final optimal joint distribution model.

[0226] The process is as follows:

[0227] 1. Review the goodness-of-fit evaluation results:

[0228] In step S4g, the log-likelihood, AIC, and BIC of each candidate Copula joint distribution model have been calculated.

[0229] 2. Preset quantitative comparison rules:

[0230] Log-likelihood value: The higher the value, the better the fit.

[0231] AIC and BIC: The smaller the value, the better the model.

[0232] AIC and BIC are usually preferred because they take into account both the complexity of the model and the complexity of the model.

[0233] 3. Compare and select the optimal model:

[0234] For each candidate model, compare its AIC and BIC values.

[0235] When AIC and BIC values ​​are similar, refer to the log-likelihood value.

[0236] The model with the smallest AIC and BIC and the highest log-likelihood value is selected as the optimal model.

[0237] Explain the relevant content of the final optimal model:

[0238] Model type: The final model selected is the R-vine model.

[0239] Model Structure: The R-vine model employs a hierarchical tree structure, containing multiple layers of Pair-Copula decomposition, which can flexibly capture complex dependencies between variables. The specific structure is as follows:

[0240] The first layer (tree T1): Using tunneling specific energy (SEC) as the core node, pair-copula connections are established with thrust (TF), cutterhead torque (TQ), and penetration (PR) respectively. Clayton Copula is used to characterize the lower tail dependency (τ). TF−SEC =0.62,τ TQ−SEC =0.58).

[0241] Second layer (tree T2): Based on the conditionalization of the first layer, establish the conditional pair-copula (C) between TF and TQ. TF,TQ∣SEC The symmetry dependency was characterized using Frank Copula (θ=3.45).

[0242] Third layer (tree T3): Establish conditional pair-copula (C) between TF and PR.TF,PR∣SEC,TQ The upper tail dependency was characterized using GumbelCopula (δ=1.78).

[0243] Goodness-of-fit metrics: The log-likelihood of the R-vine model is -1150, the AIC is 2340, and the BIC is 2443.96, indicating that it has high accuracy and low complexity when fitting data.

[0244] Parameter estimation process: The Pair-Copula parameters are solved using two-stage maximum likelihood estimation (2-stage MLE). First, the optimal marginal distributions of each variable (e.g., TF follows a Weibull distribution, parameters k=2.3, λ=2850) are transformed into uniformly distributed variables ui∈[0,1] using CDF. Then, for each Pair-Copula, the log-likelihood function is maximized. ,in Given the Copula density function, the Newton-Raphson method is used for iteration until the parameter change Δθ < 1 × 10⁻⁶. −4 .

[0245] High-energy-consuming operating conditions were identified as including:

[0246] S5a, Based on the optimal joint distribution model determined in step S4, N sets of simulated tunneling state samples are generated using the Monte Carlo simulation method. Each set of samples contains the simulated values ​​of key construction parameters and the corresponding simulated values ​​of tunneling specific energy, forming a simulated tunneling state sample set.

[0247] The process is as follows:

[0248] 1. Review the optimal joint distribution model:

[0249] In step S4, the optimal Copula joint distribution model (such as the R-vine model) has been determined.

[0250] This model describes the joint probability distribution between key construction parameters and tunneling specific energy.

[0251] 2. Monte Carlo simulation:

[0252] Sample size N: The number of samples N to be generated is determined based on statistical reliability and computational resources. Typically, N should be large enough (e.g., 10,000 or more) to ensure the stability and reliability of the sample distribution.

[0253] Generate uniformly distributed samples: Randomly generate N sets of samples from the uniform distribution U(0,1).

[0254] Inverse transformation method: Using the inverse cumulative distribution function (CDF) of the optimal joint distribution model, uniformly distributed samples are converted into actual simulated values ​​of key construction parameters and tunneling specific energy.

[0255] 3. Create a sample set of simulated tunneling conditions:

[0256] Each sample set includes simulated values ​​of key construction parameters and corresponding simulated values ​​of tunneling specific energy.

[0257] All generated samples are combined into a complete set of simulated tunneling state samples.

[0258] The 10,000 samples generated by the Monte Carlo simulation are strictly based on the optimal R-Vine model constructed in step S4. The sampling process adheres to the following constraints to ensure legality and engineering feasibility:

[0259] Data validity: All input parameters (TF, TQ, PR) are extracted from historical tunneling logs, anonymized, and do not involve personal privacy or data prohibited by law; the simulation process does not generate outliers that exceed the physical limits of the equipment (such as TF>5000kN or TQ>15MN·m).

[0260] Input-output correlation: The simulation input is a geological feature vector. With construction parameter vector The output is the tunneling specific energy. The mapping relationship strictly follows the principle of energy conservation:

[0261] .

[0262] Where D is the cutterhead diameter, v is the tunneling speed, and η is the mechanical efficiency (predicted by the model), ensuring that the output has clear physical interpretability.

[0263] Computational resource configuration: Latin hypercube sampling (LHS) is adopted to improve efficiency. The simulation time on the processor is about 45 seconds and the memory usage is <2GB, which meets the real-time requirements of the project.

[0264] S5b extracts all N simulated tunneling energy values ​​from the simulated tunneling state sample set, sorts them by numerical value, and constructs an empirical distribution sequence of simulated tunneling energy values.

[0265] The necessary procedures are as follows:

[0266] 1. Extract simulated values ​​of tunneling specific energy:

[0267] Extract the simulated tunneling specific energy (SEC) value from each sample from the set of simulated tunneling state samples generated in step S5a.

[0268] Assuming the sample set contains N sets of samples, extract N simulated values ​​of tunneling specific energy {SEC}. sim,1 SEC sim,2 ,…,SEC sim,N}

[0269] 2. Sort by numerical value:

[0270] The extracted simulated tunneling energy values ​​are sorted in ascending order to form an ordered sequence.

[0271] The sorted sequence is denoted as {SEC sorted,1 SEC sorted,2 ,…,SEC sorted,N}, where SEC sorted,1 It is the minimum value, SEC sorted,N It is the maximum value.

[0272] 3. Construct the empirical distribution sequence:

[0273] The empirical distribution sequence reflects the cumulative distribution of the simulated tunneling specific energy values.

[0274] For each sorted tunneling energy value SEC sorted,i Calculate its cumulative probability F(SEC) sorted,i )= .

[0275] The resulting empirical distribution sequence is .

[0276] S5c determines the optimal high-energy-consumption quantile based on preset statistical reliability requirements and engineering risk control requirements. The statistical reliability requirements are related to the sample size N, and the engineering risk control requirements are related to the maximum allowable risk probability of the project.

[0277] The process is as follows:

[0278] 1. Statistical reliability requirements:

[0279] Sample size N: The simulated sample size N generated in step S5a determines the accuracy of the quantile estimation. A larger N provides higher statistical reliability.

[0280] Quantile precision: Select a reasonable quantile to ensure that its estimate is statistically reliable. Generally, the selection of quantiles should avoid extreme values ​​to reduce the impact of sampling error.

[0281] 2. Project risk control requirements:

[0282] Maximum risk probability: Based on the actual needs of the project, set a maximum allowable risk probability p. max For example, if the maximum permissible risk probability for a project is 5%, then the high-energy-consumption quantile should be selected at the 95th percentile.

[0283] Quantile selection: High-energy-consuming quantile q should satisfy P(SEC>q)≤p maxThat is, the probability of exceeding this quantile does not exceed the maximum permissible risk probability.

[0284] 3. Taking into account both statistical reliability and engineering risk control:

[0285] Quantile Range: A reasonable quantile range should be selected, considering both statistical reliability and engineering risk control requirements. For example, if the sample size N=10,000, the 95% quantile has high statistical reliability while also meeting engineering risk control requirements.

[0286] Determination of the optimal quantile: Based on the above requirements, determine the optimal high-energy-consuming quantile q. best .

[0287] S5d: Based on the optimal high-energy-consumption quantile determined in step S5c, extract the corresponding high-energy-consumption judgment threshold from the empirical distribution sequence.

[0288] The necessary procedures are as follows:

[0289] 1. Review the empirical distribution sequence: In step S5b, the empirical distribution sequence of the simulated tunneling specific energy (SEC) value has been constructed, in the form of... F(SEC) sorted,i ) is the cumulative probability.

[0290] 2. Determine the optimal high-energy-consumption quantile: In step S5c, the optimal high-energy-consumption quantile q has been determined. best For example, the 95th percentile.

[0291] 3. Extract the high energy consumption threshold: In the empirical distribution sequence, find the one with the cumulative probability closest to q. best SEC (Self-Propelled Mining Energy) threshold .

[0292] Specific operation: Find the sequence of empirical distributions that satisfy F(SEC) sorted,i )≈q best SEC sorted,i .

[0293] S5e identifies samples in the simulated tunneling state sample set whose simulated tunneling energy ratio is greater than the high energy consumption judgment threshold as high energy consumption condition samples, forming a high energy consumption condition sample subset.

[0294] The necessary procedures are as follows:

[0295] 1. Review the high energy consumption threshold: In step S5d, the high energy consumption threshold SEC has been determined. threshold .

[0296] 2. Screening high-energy-consumption working condition samples: Traverse the simulated tunneling state sample set and check the simulated tunneling specific energy (SEC) value in each sample group.sim .

[0297] If the SEC sim SEC threshold If so, the sample will be identified as a high-energy-consumption operating condition sample.

[0298] 3. Form a high-energy-consumption operating condition sample subset: Collect all samples identified as high-energy-consumption operating conditions to form a high-energy-consumption operating condition sample subset.

[0299] Based on pre-set statistical reliability requirements and engineering risk control requirements, the optimal high-energy-consumption quantile is determined as follows:

[0300] S5c1, based on preset statistical reliability requirements and engineering risk control requirements, determines the search range and step size of candidate high energy consumption quantiles, and then generates a set of candidate high energy consumption quantiles.

[0301] The specific process is as follows:

[0302] 1. Determine the statistical reliability requirements:

[0303] Sample size N: The simulated sample size N generated in step S5a determines the accuracy of the quantile estimation. A larger N provides higher statistical reliability.

[0304] Quantile Range: Determine a reasonable quantile search range based on the sample size N. Generally, quantiles should avoid extreme values ​​to reduce the impact of sampling error. For example, for N=10,000, the quantile search range can be set between 90% and 99%.

[0305] 2. Determine the project risk control requirements:

[0306] Maximum risk probability p max Based on the actual needs of the project, set a maximum permissible risk probability. For example, if the maximum permissible risk probability for the project is 5%, then the high energy consumption quantile should be selected at the 95th percentile.

[0307] Quantile Range: The quantile search range is determined based on engineering risk control requirements. For example, if p max If the value is 0.05, the search range for quantiles can be set between 90% and 95%.

[0308] 3. Determine the step size:

[0309] Step size: Select an appropriate step size based on the search range of quantiles. The step size should be small enough to ensure that the optimal quantile is found accurately. For example, the step size can be set to 0.01 (i.e., 1%).

[0310] Candidate quantile set: Generate a candidate high-energy-consumption quantile set in the form {q1,q2,…,q}.m}, where q i These are candidate quantiles.

[0311] S5c2, for each candidate quantile in the candidate high-energy-consumption quantile set, perform the following operations: extract a temporary threshold from the empirical distribution sequence based on the candidate quantile, and identify a temporary high-energy-consumption sample subset from the simulated tunneling state sample set based on the temporary threshold.

[0312] The process is as follows:

[0313] 1. Review the candidate high-energy-consuming quantile set:

[0314] In step S5c1, a set of candidate high-energy-consumption quantiles has been generated, such as {0.90,0.91,0.92,0.93,0.94,0.95}.

[0315] 2. Extract temporary threshold: For each candidate quantile q i Extract the corresponding temporary threshold SEC from the empirical distribution sequence. threshold,i .

[0316] Specific steps: Find the sequence of empirical distributions whose cumulative probability is closest to q. i SEC (Self-Propelled Mining Energy) sorted,j .

[0317] 3. Identify temporary high-energy-consuming sample subsets: Use the extracted temporary threshold SEC threshold,i Samples with simulated tunneling energy values ​​greater than the threshold are selected from the simulated tunneling state sample set.

[0318] Forming a temporary high-energy-consuming sample subset S i Each sample contains simulated values ​​of key construction parameters and corresponding simulated values ​​of tunneling specific energy.

[0319] S5c3, for each temporary high-energy-consumption sample subset obtained in step S5c2, calculate its spatial clustering index in at least one preset key construction parameter dimension.

[0320] The necessary procedures are as follows:

[0321] 1. Review of the temporary high-energy-consumption sample subset:

[0322] In step S5c2, for each candidate quantile, the corresponding temporary high-energy-consuming sample subset S has been identified. i .

[0323] 2. Select key construction parameter dimensions:

[0324] Based on project requirements and data characteristics, select at least one key construction parameter dimension for analysis. For example, you can choose thrust (TF), cutterhead torque (TQ), or penetration (PR).

[0325] 3. Calculate the spatial clustering index:

[0326] For each temporary high-energy-consuming sample subset S i Calculate its spatial clustering index across selected key construction parameter dimensions. Common clustering indices include:

[0327] Mean: The average value of a sample along this dimension.

[0328] Standard deviation: The degree of dispersion of a sample along this dimension.

[0329] Coefficient of Variation (CV): The ratio of standard deviation to mean, used to measure relative dispersion.

[0330] Interquartile Range (IQR): The difference between the upper quartile and the lower quartile, reflecting the dispersion of the middle 50% of the data in a sample.

[0331] Spatial Aggregation Index (SAI): A custom metric used to measure the degree of spatial clustering of samples. For example, it can be defined by calculating the average distance between sample points.

[0332] 4. Record clustering indicators:

[0333] For each temporary high-energy-consuming sample subset S i Record its spatial clustering index in each key construction parameter dimension.

[0334] S5c4, construct a comprehensive evaluation function. The comprehensive evaluation function takes the candidate quantile as input, and its function value comprehensively reflects the statistical reliability, engineering risk control compliance and spatial clustering index calculated in step S5c3 corresponding to the candidate quantile.

[0335] The necessary procedures are as follows:

[0336] 1. Define a comprehensive evaluation function: Construct a comprehensive evaluation function F(q), where the input is the candidate high-energy-consumption quantile q, and the output is a comprehensive score reflecting the overall performance of that quantile. The function is in the form of a weighted sum:

[0337] .

[0338] Among them, w1, w2, and w3 are weighting coefficients, which are adjusted according to project requirements and importance.

[0339] 2. Calculate each component:

[0340] Statistical reliability (q): This measures the statistical reliability of a quantile q, and is typically related to the sample size N and the quantile's location. It can be measured using the standard error SE(q) of the quantile, with its reciprocal taken as the reliability metric.

[0341] .

[0342] RiskControl(q): Measures whether the quantile q meets the maximum allowable risk probability p of the project. max Calculate the probability P(SEC>q) corresponding to the quantile and evaluate whether it satisfies p. max :

[0343] .

[0344] Spatial aggregation index Aggregation(q):

[0345] Use a spatial clustering index calculated in step S5c3, such as the coefficient of variation (CV) or the spatial clustering index (SAI). For example, use the reciprocal of the coefficient of variation as a clustering index:

[0346] .

[0347] 3. Weighting:

[0348] Assign weights w1, w2, and w3 based on project requirements and importance. For example:

[0349] w1=0.4 (statistical reliability).

[0350] w2=0.3 (Engineering risk control compliance).

[0351] w3=0.3 (spatial clustering index).

[0352] 4. Calculate the comprehensive evaluation value: For each candidate quantile q, calculate its comprehensive evaluation value F(q):

[0353] .

[0354] S5c5 substitutes each quantile in the candidate high-energy-consumption quantile set into the comprehensive evaluation function for calculation, and selects the candidate quantile that makes the function value optimal as the best high-energy-consumption quantile.

[0355] The necessary procedures are as follows:

[0356] 1. Review the candidate high-energy-consuming quantile set:

[0357] In step S5c1, a set of candidate high-energy-consumption quantiles has been generated, such as {0.90,0.91,0.92,0.93,0.94,0.95}.

[0358] 2. Calculate the overall assessment value for each quantile:

[0359] For each candidate quantile q i The comprehensive evaluation value F(q) is calculated using the comprehensive evaluation function F(q) constructed in step S5c4. i ).

[0360] 3. Compare the overall evaluation values:

[0361] Compare the overall evaluation values ​​F(q) of all candidate quantiles i We select the quantile that yields the optimal function value. Typically, the optimal quantile corresponds to the maximum comprehensive evaluation value.

[0362] 4. Determine the optimal high-energy-consuming quantile:

[0363] The candidate quantile with the highest comprehensive evaluation value was selected as the optimal high-energy-consuming quantile q. best .

[0364] Calculate the spatial clustering index of the temporary high-energy-consumption sample subset on at least one preset key construction parameter dimension, including:

[0365] S5c3a. Extract at least one pair of preset key construction parameters from the temporary high-energy-consumption sample subset to form a sample point set in a two-dimensional parameter space.

[0366] The necessary procedures are as follows:

[0367] 1. Review of the temporary high-energy-consumption sample subset:

[0368] In step S5c2, for each candidate high-energy-consumption quantile, a corresponding temporary high-energy-consumption sample subset S has been identified. i Each sample subset contains simulated values ​​of key construction parameters and tunneling specific energy.

[0369] 2. Select key construction parameter pairs:

[0370] Based on project requirements and data characteristics, select at least one pair of key construction parameters. These parameter pairs are generally considered to have a significant impact on energy consumption in engineering practice. Common parameter pairs include:

[0371] Thrust (TF) and cutterhead torque (TQ);

[0372] Thrust (TF) and penetration (PR);

[0373] Cutter torque (TQ) and penetration (PR);

[0374] 3. Extracting the Sample Point Set: Extract the values ​​of selected key construction parameter pairs from the temporary high-energy-consumption sample subset to form a sample point set in a two-dimensional parameter space. Each sample point is a two-dimensional vector representing the value of a pair of key construction parameters.

[0375] S5c3b: The kernel density estimation algorithm is used to estimate the sample point set to obtain the joint probability density distribution of the two-dimensional parameter space.

[0376] The necessary procedures are as follows:

[0377] 1. Review the sample point set: In step S5c3a, the sample point set in the two-dimensional parameter space has been extracted. For example, if the thrust force (TF) and cutterhead torque (TQ) are selected as the key construction parameter pair, the sample point set is: {(2800,9.2),(3000,10.5),(2900,9.8),…}.

[0378] 2. Choosing a kernel density estimation algorithm: Kernel density estimation is a nonparametric estimation method used to estimate the probability density function of data. Commonly used kernel functions include Gaussian kernel, uniform kernel, and triangular kernel. The Gaussian kernel is the most commonly used choice because it has good smoothness and mathematical properties.

[0379] 3. Select bandwidth parameters:

[0380] The bandwidth parameter (h) is an important parameter in kernel density estimation, controlling the smoothness of the kernel function. Methods for selecting the bandwidth include:

[0381] Silverman's rule: Applicable to the normal reference method, the calculation formula is:

[0382] , where n is the sample size and std(X) is the standard deviation of the sample.

[0383] Scott's rule: Applicable to large sample sizes, the calculation formula is:

[0384] .

[0385] 4. Calculate the joint probability density distribution:

[0386] For each sample point (x) i ,y i The joint probability density is calculated using the kernel function K and bandwidth h:

[0387] .

[0388] Where K is the kernel function, usually a Gaussian kernel is chosen:

[0389] .

[0390] 5. Generate joint probability density map: Generate a grid in the two-dimensional parameter space and calculate the joint probability density value at each grid point.

[0391] Visualize the joint probability density distribution and generate a two-dimensional probability density map.

[0392] S5c3c: Calculate the peak density of the joint probability density distribution and use it as a spatial clustering index to characterize the spatial clustering degree of temporary high-energy-consuming sample subsets.

[0393] The process is as follows:

[0394] 1. Review the joint probability density distribution:

[0395] In step S5c3b, the joint probability density distribution of the two-dimensional parameter space has been obtained through the kernel density estimation algorithm. .

[0396] 2. Determine the peak density:

[0397] In the joint probability density distribution, find the maximum value, i.e., the peak density f. peak Peak density is the highest point in the joint probability density distribution, representing the region where the sample points are most concentrated.

[0398] 3. Calculate peak density:

[0399] Traverse all grid points of the joint probability density distribution and find the maximum value:

[0400] f peak =max{ Grid points}

[0401] 4. Record the peak density: Record the peak density f peak Record it as an indicator of the spatial clustering of this temporary high-energy-consuming sample subset.

[0402] Step S6 includes:

[0403] S6a. Extract at least one pair of key construction parameters from the high-energy-consumption working condition samples to form a two-dimensional risk sample point set.

[0404] The necessary procedures are as follows:

[0405] 1. Review of high-energy-consumption working condition samples: In step S5e, a subset of high-energy-consumption working condition samples has been identified. These samples represent the high-energy-consumption state during shield tunneling.

[0406] 2. Select a pair of critical construction parameters: Based on project requirements and data characteristics, select at least one pair of critical construction parameters. These parameter pairs are generally considered to have a significant impact on energy consumption and risk in engineering practice. Common parameter pairs include:

[0407] Thrust (TF) and cutterhead torque (TQ);

[0408] Thrust (TF) and penetration (PR);

[0409] Cutter torque (TQ) and penetration (PR);

[0410] 3. Extract the two-dimensional risk sample point set:

[0411] The values ​​of selected key construction parameter pairs are extracted from the high-energy-consumption operating condition sample subset to form a two-dimensional risk sample point set. Each sample point is a two-dimensional vector representing the value of a pair of key construction parameters.

[0412] S6b. The kernel density estimation algorithm is used to estimate the two-dimensional risk sample point set to obtain the two-dimensional joint probability density distribution, and a probability density heatmap is generated based on the distribution.

[0413] The necessary procedures are as follows:

[0414] 1. Review the two-dimensional risk sample point set: In step S6a, the two-dimensional risk sample point set has been extracted. For example, if the thrust force (TF) and cutterhead torque (TQ) are selected as the key construction parameter pair, the sample point set is: {(2800,9.2),(3000,10.5),(2900,9.8),…}.

[0415] 2. Select the kernel density estimation algorithm:

[0416] Kernel density estimation is a nonparametric estimation method used to estimate the probability density function of data. Commonly used kernel functions include Gaussian kernel, uniform kernel, and triangular kernel. The Gaussian kernel is the most commonly used choice because it has good smoothness and mathematical properties.

[0417] 3. Select bandwidth parameters:

[0418] The bandwidth parameter (h) is an important parameter in kernel density estimation, controlling the smoothness of the kernel function. Methods for selecting the bandwidth include:

[0419] Silverman's rule: Applicable to the normal reference method, the calculation formula is:

[0420] Where n is the sample size and std(X) is the standard deviation of the sample.

[0421] Scott's rule: Applicable to large sample sizes, the calculation formula is:

[0422] .

[0423] 4. Calculate the two-dimensional joint probability density distribution:

[0424] For each sample point (x) i ,y i The joint probability density is calculated using the kernel function K and bandwidth h:

[0425] .

[0426] Where K is the kernel function, usually a Gaussian kernel is chosen: .

[0427] h x and h y These are the bandwidths in the x and y directions, respectively.

[0428] 5. Generate probability density heatmap: Generate a grid in the two-dimensional parameter space and calculate the joint probability density value at each grid point.

[0429] Visualize the joint probability density distribution and generate a two-dimensional probability density heatmap. Darker colors in the heatmap indicate higher probability density, meaning the sample points are more concentrated.

[0430] S6c: Based on the two-dimensional joint probability density distribution, extract probability density contour lines, and determine multiple candidate high-risk areas according to the preset contour line density value threshold.

[0431] The necessary procedures are as follows:

[0432] 1. Review of the two-dimensional joint probability density distribution: In step S6b, the joint probability density distribution of the two-dimensional parameter space has been obtained through the kernel density estimation algorithm. .

[0433] 2. Determine contour density thresholds: Based on engineering requirements and risk assessment standards, preset one or more contour density thresholds. These thresholds are used to distinguish between high-risk and low-risk areas. For example, the 90th, 95th, and 99th percentiles of the probability density distribution can be selected as thresholds.

[0434] 3. Extracting Probability Density Contours: In a two-dimensional parameter space, extract contour lines whose probability density equals a preset threshold. A contour line is the set of all points with the same density value in a probability density distribution. Use a contour extraction algorithm (such as the `contour` function in matplotlib) to generate a contour plot.

[0435] 4. Identify candidate high-risk areas: Based on a preset contour density threshold, identify areas with a probability density higher than this threshold. These areas correspond to high-risk parameter combinations. These areas can be visually identified using a contour map, and their extent can be recorded.

[0436] S6d. Spatial clustering and merging of multiple candidate high-risk regions are performed to form the final high-risk parameter combination interval. The spatial clustering and merging process is based on the spatial distance and density similarity between regions.

[0437] The necessary procedures are as follows:

[0438] 1. Review of candidate high-risk areas:

[0439] In step S6c, multiple candidate high-risk regions have been identified, each defined by a range in the parameter space. For example:

[0440] High-risk area 1: 2900≤TF≤3000, 9.5≤TQ≤10.5.

[0441] High-risk area 2: 2950≤TF≤3050, 10.0≤TQ≤11.0.

[0442] 2. Define spatial distance and density similarity:

[0443] Spatial distance: Calculates the Euclidean distance between the center points of two regions.

[0444] Density similarity: Compare the probability density distributions of two regions to ensure that they are similar in density.

[0445] 3. Spatial clustering:

[0446] Clustering algorithms (such as DBSCAN and K-Means) are used to cluster candidate high-risk regions. DBSCAN is a density-based clustering algorithm that is suitable for identifying high-density regions with different shapes and sizes.

[0447] For each region, its center point and density features are calculated, and then adjacent and similar regions are grouped into the same cluster based on spatial distance and density similarity.

[0448] 4. Merging process:

[0449] For regions within the same cluster, a merging process is performed. During merging, the union of the parameter ranges of all regions is taken to form a larger interval of high-risk parameter combinations.

[0450] For example, if two regions overlap or are adjacent in TF and TQ, their ranges are merged.

[0451] S6e. Compare the final high-risk parameter combination range with the preset construction safety parameter range, eliminate the range that exceeds the feasible range of the project, and output the high-risk parameter combination range verified by the project for construction risk warning and operation optimization.

[0452] The necessary procedures are as follows:

[0453] 1. Review the final high-risk parameter combination interval: In step S6d, the final high-risk parameter combination interval was formed through spatial clustering and merging. For example:

[0454] High-risk parameter combination range: 2900≤TF≤3050, 9.5≤TQ≤11.0.

[0455] 2. Preset construction safety parameter range:

[0456] Based on engineering design and safety standards, preset safe ranges for construction parameters. For example:

[0457] Thrust (TF) safety range: 2500≤TF≤3200.

[0458] Cutter head torque (TQ) safety range: 8.0≤TQ≤12.0.

[0459] 3. Comparison and verification:

[0460] Check whether the high-risk parameter combination range is completely within the preset construction safety parameter range.

[0461] If the parameter range of a certain interval exceeds the preset safe range, then the interval will be removed or its range will be adjusted to the safe range.

[0462] 4. Output the high-risk parameter combination range after engineering verification:

[0463] The high-risk parameter combination range after verification is used as the final result for construction risk early warning and operation optimization.

[0464] The embodiments described in this specific implementation are preferred embodiments of this application and are not intended to limit the scope of protection of this application. Therefore, all equivalent changes made in accordance with the structure, shape and principle of this application should be covered within the scope of protection of this application.

Claims

1. A method for energy consumption modeling and multi-parameter risk identification in shield tunneling construction in complex geological formations, characterized in that, include: S1: Collect construction parameters and corresponding geological description information during the shield tunneling process, and extract structured features from the geological description information to generate geological feature vectors; S2. Calculate the correlation between each construction parameter and the preset energy consumption index, and based on the correlation analysis results, select key construction parameters from the construction parameters; S3. Using geological feature vectors and key construction parameters as inputs, and tunneling specific energy as the prediction target, construct a Gaussian process regression model, and optimize the model by combining multiple preset kernel functions to establish an energy consumption prediction model; S4. Based on the key variables constituting the input and output of the energy consumption prediction model, fit their marginal probability distributions respectively, and on this basis, construct Copula joint distribution models with different dependency structures for comparison, and then select the optimal joint distribution model according to the goodness-of-fit criterion; S5. Use the optimal joint distribution model to perform Monte Carlo simulation, generate simulated tunneling state samples, and set a high energy consumption threshold according to the statistical distribution of tunneling specific energy in the samples to identify high energy consumption working condition samples; S6. Extract key construction parameter combinations from the high energy consumption working condition samples, calculate their joint probability density and generate a probability density heatmap to identify high-risk parameter combination intervals for construction risk warning and operation optimization.

2. The method for energy consumption modeling and multi-parameter risk identification in shield tunneling construction in complex strata according to claim 1, characterized in that, Using geological feature vectors and key construction parameters as inputs, and tunneling specific energy as the prediction target, a Gaussian process regression model is constructed. This model is then optimized by combining various pre-defined kernel functions to establish an energy consumption prediction model. The process includes: S31, selecting multiple basic kernel functions from a pre-defined kernel function set, including at least radial basis function kernels, Matern kernels, rational quadratic kernels, and white noise kernels; S32, combining the basic kernel functions in pairs to generate multiple hybrid kernel function configurations; S33, using geological feature vectors and key construction parameters as inputs, and tunneling specific energy as the prediction target, constructing and training corresponding Gaussian process regression sub-models using each hybrid kernel function configuration; S34, evaluating the prediction performance of each sub-model based on a validation dataset, with evaluation metrics including the coefficient of determination and root mean square error; and S35, selecting the Gaussian process regression sub-model with the best prediction performance from all sub-models as the final energy consumption prediction model.

3. The method for energy consumption modeling and multi-parameter risk identification in shield tunneling construction in complex strata according to claim 2, characterized in that, The process of combining basic kernel functions in pairs to generate multiple hybrid kernel function configurations includes: S321, generating all possible unordered pairwise combinations based on the selected basic kernel functions; S322, generating a corresponding hybrid kernel function configuration for each unordered pairwise combination through kernel function addition or kernel function multiplication.

4. The method for energy consumption modeling and multi-parameter risk identification in shield tunneling construction in complex strata according to claim 1, characterized in that, Based on the key variables constituting the input and output of the energy consumption prediction model, the marginal probability distributions are fitted respectively, including: S4a, determining the set of key variables to be fitted to the marginal probability distribution, which includes the geological feature vector and key construction parameters constituting the input of the energy consumption prediction model, and the tunneling specific energy constituting its output; S4b, for each variable in the set of key variables, selecting multiple candidate distributions from a preset candidate probability distribution family for fitting, the candidate probability distribution family includes multiple types of normal distribution, log-normal distribution, Weibull distribution, Beta distribution, Logistic distribution, and Gumbel distribution; S4c, evaluating the fitting effect of each candidate distribution fitted for each variable in step S4b based on at least one goodness-of-fit index among the Akaike information criterion, Bayesian information criterion, and log-likelihood value; S4d, for each variable, selecting the optimal marginal probability distribution from the multiple candidate distributions fitted to it according to the goodness-of-fit index.

5. The method for energy consumption modeling and multi-parameter risk identification in shield tunneling construction in complex strata according to claim 4, characterized in that, For each variable in the set of key variables, multiple candidate distributions are selected from a pre-defined family of candidate probability distributions for fitting, including: S4b1, based on the data distribution characteristics of each variable in the set of key variables, one or more corresponding candidate probability distribution types are dynamically matched or selected from the pre-defined family of candidate probability distributions for each variable; the data distribution characteristics are determined at least by skewness, kurtosis, and empirical distribution shape; S4b2, for each candidate probability distribution type selected for each variable, the parameters are estimated using the maximum likelihood estimation method or the method of moments based on the historical observation data of that variable, and the fitting of that candidate distribution is completed.

6. The method for energy consumption modeling and multi-parameter risk identification in shield tunneling construction in complex strata according to claim 5, characterized in that, Based on this, Copula joint distribution models with different dependency structures are constructed and compared. Then, the optimal joint distribution model is selected according to the goodness-of-fit criterion, including: S4e, constructing at least two candidate Copula joint distribution models with different dependency structures in parallel, wherein the dependency structure type is selected from a set consisting of C-vine, D-vine, and R-vine structures; S4f, for each candidate Copula joint distribution model constructed in step S4e, based on the optimal marginal probability distribution selected for each key variable in step S4d, using the maximum likelihood estimation method to solve for the P-value of the model. The air-Copula function parameters are used to complete the model parameter estimation; S4g, based on a unified set of goodness-of-fit criteria, quantifies and evaluates the overall goodness-of-fit of each candidate Copula joint distribution model whose parameters have been estimated in step S4f. The set of goodness-of-fit criteria includes the Akaike information criterion, the Bayesian information criterion, and the log-likelihood value; S4h, based on a preset quantitative comparison rule, compares the evaluation results of each candidate Copula joint distribution model according to the set of goodness-of-fit criteria obtained in step S4g, selects the candidate model with the best comprehensive evaluation result, and determines it as the final optimal joint distribution model.

7. The method for energy consumption modeling and multi-parameter risk identification in shield tunneling construction in complex strata according to claim 1, characterized in that, The identification of high-energy-consuming working condition samples includes: S5a, based on the optimal joint distribution model determined in step S4, generating N sets of simulated tunneling state samples using the Monte Carlo simulation method. Each set of samples contains simulated values ​​of key construction parameters and corresponding simulated values ​​of tunneling specific energy, forming a simulated tunneling state sample set; S5b, extracting all N simulated values ​​of tunneling specific energy from the simulated tunneling state sample set, sorting them by numerical value, and constructing an empirical distribution sequence of the simulated values ​​of tunneling specific energy; S5c, determining the optimal high-energy-consuming quantile based on preset statistical reliability requirements and engineering risk control requirements, where the statistical reliability requirement is related to the sample size N, and the engineering risk control requirement is related to the maximum allowable risk probability of the project; S5d, extracting the corresponding high-energy-consuming judgment threshold from the empirical distribution sequence according to the optimal high-energy-consuming quantile determined in step S5c; S5e, identifying samples in the simulated tunneling state sample set whose simulated values ​​of tunneling specific energy are greater than the high-energy-consuming judgment threshold as high-energy-consuming working condition samples, forming a high-energy-consuming working condition sample subset.

8. The method for energy consumption modeling and multi-parameter risk identification in shield tunneling construction in complex strata according to claim 7, characterized in that, Based on preset statistical reliability requirements and engineering risk control requirements, the optimal high-energy-consumption quantile is determined as follows: S5c1, based on preset statistical reliability requirements and engineering risk control requirements, the search range and step size of candidate high-energy-consumption quantiles are determined, thereby generating a set of candidate high-energy-consumption quantiles; S5c2, for each candidate quantile in the set of candidate high-energy-consumption quantiles, the following operations are performed: a temporary threshold is extracted from the empirical distribution sequence based on the candidate quantile, and a temporary high-energy-consumption sample subset is identified from the simulated tunneling state sample set based on the temporary threshold; S5c3, for step S... For each temporary high-energy-consumption sample subset obtained in step 5c2, calculate its spatial clustering index in at least one preset key construction parameter dimension; in step S5c4, construct a comprehensive evaluation function, which takes candidate quantiles as input, and its function value comprehensively reflects the statistical reliability, engineering risk control compliance, and spatial clustering index calculated in step S5c3 corresponding to the candidate quantile; in step S5c5, substitute each quantile in the candidate high-energy-consumption quantile set into the comprehensive evaluation function for calculation, and select the candidate quantile that makes the function value optimal as the best high-energy-consumption quantile.

9. The method for energy consumption modeling and multi-parameter risk identification in shield tunneling construction in complex strata according to claim 8, characterized in that, Calculating the spatial clustering index of a temporary high-energy-consumption sample subset in at least one preset key construction parameter dimension includes: S5c3a, extracting at least one pair of preset key construction parameters from the temporary high-energy-consumption sample subset to form a sample point set in a two-dimensional parameter space; S5c3b, estimating the sample point set using a kernel density estimation algorithm to obtain the joint probability density distribution of the two-dimensional parameter space; S5c3c, calculating the peak density of the joint probability density distribution and using it as a spatial clustering index characterizing the degree of spatial clustering of the temporary high-energy-consumption sample subset.

10. The method for energy consumption modeling and multi-parameter risk identification in shield tunneling construction in complex strata according to claim 1, characterized in that, Step S6 includes: S6a, extracting at least one pair of key construction parameters from the high-energy-consumption working condition sample to form a two-dimensional risk sample point set; S6b, estimating the two-dimensional risk sample point set using a kernel density estimation algorithm to obtain a two-dimensional joint probability density distribution, and generating a probability density heatmap based on this distribution; S6c, extracting probability density contour lines based on the two-dimensional joint probability density distribution, and determining multiple candidate high-risk areas according to a preset contour line density value threshold; S6d, performing spatial clustering and merging processing on the multiple candidate high-risk areas to form the final high-risk parameter combination interval, wherein the spatial clustering and merging processing is based on the spatial distance and density similarity between areas; S6e, comparing the final high-risk parameter combination interval with the preset construction safety parameter range, eliminating intervals that exceed the feasible range of the project, and outputting the high-risk parameter combination interval verified by the project for construction risk early warning and operation optimization.

Images