Frontiers in Computing and Intelligent Systems ISSN: 2832-6024 | Vol. 8, No. 3, 2024 6 Analysis of Natural Gas Hydrate Resource Content based on Kriging Algorithm and TOPSIS Evaluation Model Boyang Guan * Mechanical Engineering College, Xi’an Shiyou University, Xi’an, Shaanxi, China * Corresponding author: Boyang Guan Abstract: The exploration and evaluation system for natural gas hydrate resources is currently incomplete. Research on evaluation methods for the quantity of natural gas hydrate resources is an important part of theoretical studies and is crucial in guiding exploration and development activities. In this paper, we establish the total effective thickness (Z), average porosity of the reservoir ( ), and the absolute value of average saturation ( ) as indicators for evaluating the quantity of natural gas hydrate resources. We then statistically analyze the data from various exploration points to obtain the values of these indicators. Subsequently, we employ the Kriging algorithm to interpolate the data, obtaining a more densely distributed grid. Following this, we establish a TOPSIS evaluation model to calculate the scores of each exploration point. Finally, a heat map is generated to represent the distribution of resource quantity scores in various regions. The results indicate that natural gas hydrate resources are primarily distributed in the upper left, lower right, and directly below the detection area. Keywords: Natural Gas Hydrate; Kriging Algorithm; TOPSIS Evaluation Model. 1. Introduction At present, the exploration and evaluation system of natural gas hydrate resources is still imperfect, and the research on the evaluation method of natural gas hydrate resources is an important part of the theoretical research on natural gas hydrate, which is the key to effectively guide the exploration and development activities of natural gas hydrate. At present, some scholars have already carried out related research. xiao-Han Liu et al [1] used the newly proposed genesis method and geological analogy method to evaluate the gas hydrate resources. The results show that the conventional oil and gas resources in the South China Sea are 346.29 × 108 t, and the volumes of gas hydrate and free power field are 25.19 × 104 km3 and (2.05 ~ 2.48) × 106 km3, respectively, and the total in-situ gas hydrate resources are about (4.47 ~ 6.02) × 1012 m3. Xing-Wen Zhang et al [2] divided the South China Sea natural gas hydrate investigation into seven stages, based on previous studies, obtaining the key geological parameters of each stage by analogy with other regions, evaluating the natural gas hydrate resources in these seven stages by using the volumetric method, and simulating the downward trend of these estimates by using a new trend analysis method, and finally predicting the natural gas hydrate resources in the South China Sea in 2025 and 2030. Qing - guo Meng et al [3] conducted hydrate formation experiments using core samples of hydrate-containing sediments from the Qilian Mountain permafrost zone and multicomponent natural gas with similar composition to the natural gas hydrate from the Qilian Mountain permafrost zone. The formation and distribution characteristics of multicomponent natural gas hydrate in the core samples under high-pressure and low- temperature conditions were observed in situ using X-ray computed tomography (X-CT). Enze Wang et al [4] established a three-dimensional evaluation model for the quality of shale oil resources based on multi-stage pyrolysis parameters. In order to improve the applicability of the model, a simple total oil yield (TOY) algorithm was proposed, which can be used for rapid and accurate evaluation of shale resource quality where only conventional pyrolysis data are available. Jinbu Li et al [5] proposed a fast method to evaluate the total oil yield using a single conventional rock evaluation experiment. Firstly, temperature thresholds (TOK) for s2 oil and cracked hydrocarbons were determined by combining a char plot of the received shale with solvent-extracted replicates. Then, total yields were obtained directly from hydrocarbon evaporation at temperatures below TOK in conventional rock experiments. The results show that the TOK value is controlled by the maturity and pore structure of the sample. The higher the maturity, the higher the specific surface area, the smaller the pore size and the higher the TOK. The exploration and evaluation system for natural gas hydrate resources remains incomplete, necessitating the development of effective methods for evaluating their quantity. This research constitutes a vital component of theoretical studies on natural gas hydrates and is essential in guiding exploration and development activities. In this study, we address this need by establishing indicators for resource quantity evaluation and employing the Kriging algorithm for data interpolation to achieve a more comprehensive dataset. Subsequently, a TOPSIS evaluation model is established to compute scores for each exploration point, providing insights into the distribution of natural gas hydrate resources within various regions. In this paper, using http://www.nmmcm. org.cn/uploads/20240510/2024%E5%B9%B4%E7%AC%A C%E4%B9%9D%E5%B1%8A%E6%95%B0%E7%BB%B 4%E6%9D%AF%E6%95%B0%E5% AD%A6%E5% BB% BA%E6%A8%A1%E6%8C%91%E6%88%98%E8%B5%9 B%E9%A2%98%E7%9B%AE.rar to describe the algorithms used in this paper. 2. Construction and Calculation of Evaluation Indicators In order to evaluate the gas hydrate resource of an 7 exploration site, we constructed three key metrics to characterize it, including total effective thickness ( ), mean porosity ( ) and mean absolute saturation ( ). These metrics will be used as the basis for assessing the gas hydrate resource potential of this exploration site. These indicators are described in detail below: 1) Total effective thickness ( ) In geological exploration, total effective thickness is a key geological parameter that describes the sum of the vertical thicknesses of effective rock layers with gas hydrate potential in the stratigraphy in which the exploration site is located. We have defined this metric in order to quantify the distribution of rock layers available for gas hydrate reservoirs. In the dataset given in this question, the thickness was counted for each successive depth series with non-zero saturation and summed cumulatively. 2) Average formation porosity ( ) The average porosity of an oil formation is the proportion of the total volume occupied by the pore space in the rock or sediment in the formation at the point where the oil formation occurs. The level of average porosity will affect the permeability and storage capacity of the reservoir and is therefore crucial for assessing the potential storage capacity of the exploration site. In the dataset given in this question, the average value of porosity is calculated for the non-0 saturation points of each exploration site. 3) Average saturation absolute value ( ) The mean saturation absolute value is the average value of the absolute saturation of gas hydrate in a gas hydrate reservoir. This index is calculated by a combination of seismic data, core analysis and numerical simulation. The absolute saturation value reflects the degree of saturation of gas hydrate in the reservoir and provides an important quantitative description of the resource potential of the exploration site. In the dataset given in this question, the average value is calculated after taking the absolute value of non-zero saturation for each exploration point. The calculated metrics for each exploration point are shown in Table 1. Table 1. Characteristic values for each exploration site exploration site Total Effective Thickness Average Porosity of Oil Formation The absolute value of average saturation W01 17.8308 0.497964957 0.017823932 W02 136.2 0.520521855 0.183099285 W03 30.4 0.457147039 0.018322039 W04 101.4 0.504006011 0.261930773 W05 77.5 0.48971977 0.055410654 W06 135.4 0.465867998 0.043989583 W07 78.6384 0.533092054 0.102083915 W08 126.492 0.484948795 0.214307669 W09 154.2288 0.501111759 0.197240909 W10 55.1688 0.527341989 0.135443094 W11 103.9 0.427846867 0.042548427 W12 16.1544 0.536748113 0.007519811 W13 134 0.542675896 0.05914505 W14 97.536 0.487260156 0.020013281 3. Spatial Interpolation based on Kriging Interpolation Algorithm 3.1. Modelling Krige's interpolation, also known as spatial self-covariance optimal interpolation, is a method of linearly optimal, unbiased interpolation estimation of spatially distributed data. The method was by the South African engineer Krige (Krige D G) and statistician Sichel (Sichel H S) in the 1950s according to the samples of different spatial locations and different degrees of correlation between the samples, each sample is given a certain weight, a sliding weighted average, to estimate the average value of the samples on the unknown sample point of a method, proposed by the French statistician Matheron. It commemorates the first application of statistical techniques to geomorphological assessment by Kriging in 1951. It is a method of linear unbiased optimal estimation of the value of a regionalised variable at an unsampled point, using the raw data of the regionalised variable and the structural characteristics of the variogram. This method is an unbiased optimal estimation of the sample points to be estimated on the basis of analysing the shape, size, spatial location and interrelationships of the measured sample points, the mutual spatial location of the measured sample points and the sample points to be estimated, and the structural information provided by the variogram function. From the interpolation point of view, this method is applicable under the condition that there is a spatial correlation of the regionalized variables, and it makes maximum use of the information provided by spatial sampling, not only considering the data of the sample point, but also considering the data of the neighbouring sample points, not only considering the spatial position of the sample point to be estimated and the neighbouring known sample points, but also considering the positional relationship between the neighbouring sample points and the structural characteristics of the spatial distribution of the existing observations, thus making the estimation of the sample point to be estimated more optimum and efficient. The structural characteristics of the distribution of the existing observations, so that this interpolation method is more accurate than other methods, and can give the estimation error. Kriging (Kriging) interpolation method contains three important concepts, one is the regionalisation of variables; the second is the variance function; the third is the covariance function. 1) Theory of regionalised variables A variable is spatially related to its location, i.e. when a variable is spatially distributed, the variable is called a regionalised variable. This kind of variable reflects the distribution characteristics of a certain attribute in space, and the natural gas hydrate resource problem solved in this question also has clear spatial attributes. The regionalised variable has a dual nature; before observation, the regionalised variable Z(X) is a random field and after observation, it is a definite spatial point function value. The number of interpolated points and the density of interpolated points are examples of regionalised variables, and their spatial structure and spatial dependence relationships can be investigated using regionalised variable theory. 2) Variance function The variational function, also known as the variational function and moments of variation, studies the spatial correlation and spatial structure of variables by determining the difference between regionalised variables separating equidistant sample points. In one-dimensional conditions the variational function is defined as follows: when the spatial point x varies on the one-dimensional x-axis, the values of the regionalization variable Z(x) at the points x and x+h are Z(x) and Z(x+h), and the regionalization variable Z(x) has a finite 8 variance in the increment for all the distance vectors h, and it does not depend on the position x. Therefore, the regionalization variable Z(x) in the separation distance h between two points x and x+h can be expressed in terms of the mathematical expectation of the squares of their increments (i.e., the variance of the increments of the regionalised variables) as follows. 1 2 var (1) r(h) is called the semi-variational function. The semi- variogram is a plot of the semi-variogram r(h) value as a function of distance h. It has three characteristic parameters: the abutment value, the range, or spatial dependence range, and the nugget value, or regional discontinuity value, which determine the properties of the semi-variogram. For the observed data series z(x) (where i = 1,2,3,... n ), the value of the sample semi-variance function can be calculated by the following equation. 1 2 ∑ (2) where: N(h) is the logarithm of the data pair , , , split by h; z(x,) and z(x+h) are the measurements of the sample at points x, and x+h, respectively; h is the distance between sample points. The variability function reveals the pattern of spatial variability over the entire scale, and the variability function is only meaningful at the maximum interval distance of 1/2. 3) Covariance function The Kriging method is a geostatistical method for estimating the interpolated values between observation samples, which is based on the theory of regionalised variables, and when a simulation model of the semi-variance function of a given variable has been obtained, the sample point observations can be used to estimate the value of the regionalised variable at an unsampled point on the study area. Let the regionalised variable be z(x), let z(x0) be the estimated value of any point x, near which there are n known observation samples x,(i= 1,2,3,...,n), and zx,) be the observed value of the variable at the sample point x, if it is assumed that the regionalised variable z(x) satisfies the second-order smoothness assumption and the intrinsic assumption and that it has the mathematical expectation of m and the covariance of c(h). A variation function of r(h) exists, i.e... 1 2 (3) Then the estimate at the prediction point x0. can be found by a linear combination of the observations at n surrounding sampling points, i.e.: (4) The formula is the weight of the sampling point , and the weight of the distance inverse interpolation method is different from the weight of , the assigned weight value not only takes into account the distance between the prediction point and the sampling point but also takes into account the spatial distribution of the prediction point and the sampling point and the sampling point between each other. To obtain the unbiased optimal estimate, the following two conditions must be satisfied. (1) Unbiased estimation, i.e. ′ 0 (2) Minimal estimation variance, i.e. ′ 一 → Conditional on the existence of a variational function, the ordinary Kriging equation and the variance of Kriging estimates are expressed in terms of the variational function according to the relation between the covariance and the variational function: as.   , , , 1, 1,2,⋯ (5) , , (6) and the weights, , are required to satisfy the above equations. Where the equation , is the semi-variance between the observation point x, and x,, , is the semi-variance between the sampling point x_i and the interpolated point x, and μ is the Lagrange multiplier related to the minimisation of variance. From this equation, the values of the weights are calculated, and the interpolated value z(x) at the point x to be estimated can be found by substituting into the above equation. The equation can also be expressed as a matrix as: 9 ⋯ 1 ⋯ 1 ⋮ ⋮ ⋮ ⋮ ⋮ ⋯ 1 1 1 ⋯ 1 0 , 2 ⋮ , , , ⋮ , 1 (7) Then the Kriging equation is: (8) Its estimated variance is: , (9) 3.2. Solution of the Model The idea of Kriging interpolation is that it first considers the distribution of variability of spatial attributes in spatial locations, determines the range of distances that have an effect on the value of a to-be-interpolated point, and then estimates the value of the to-be-interpolated point attribute by using the sampling points within this range, and it is a method for finding the optimal linear and unbiased interpolation estimator. It is a method of estimating the block grade by assigning certain coefficients to each sample value separately and finally performing a weighted average in order to achieve linear, unbiased and minimum estimation variance estimation after considering the geometric features such as the shape and size of the information samples and their to-be-estimated block segment's spatial distribution position among each other as well as the spatial structure of the taste. Its calculation steps are: 1) Input raw data, i.e. sampling points. 2) Gridding, selecting the extent of the area and the size of the grid. 3) Data inspection and analysis, observation to see whether the sampling value is in line with the actual situation, deleting the obvious difference points. 4) Calculation of histograms in order to decide whether to pre-process the raw data. 5) Calculation of the variational function using the variational function principle to understand the spatial structure of the variables. 6) Kriging interpolation estimation 1) Estimation of weight coefficients of points to be estimated Using the method of polygonal estimation, first determine the weight of the sampling point closest to the point to be estimated, and then carry out the estimation of the weight of the sampling point according to the following formula. 1 ∑ (10) After getting the weights of the nearest sampling points, determine the weights of the points to be estimated according to the method of determining the weights of the points to be estimated. 2) According to the search strategy, select the appropriate reference estimation point. 3) According to the variance function that has been derived and the number of sampling points, assuming that there are three sampling points, three equations are listed to find the coefficients of the system of equations, which is given by the formula: 1,1 1,2 1,3 2,1 2,2 2,3 3,1 3,2 3,3 0,1 0,2 0,3 (11) 4) Analyse the effect on the weight values of changing the nugget value under the condition of isotropy versus changing the anisotropy under the condition of the same nugget value. 5) According to the derived weight value, it can be substituted into the above equation to find the linear combination of n sampling values in the assessment domain. The interpolation results are shown in Fig. 1, Fig. 2 and Fig. 3. Fig 1. Total Effective Thickness Interpolation Results 10 Fig 2. Interpolation results of the average porosity of the oil formation Fig 3. Average saturation absolute value interpolation results 4. Evaluation of Gas Hydrate Resource Content at Exploration Well Locations based on TOPSIS Algorithm From the equation of the linear relationship between resources and reservoir parameters, it can be seen that gas hydrate resources are affected by effective area, effective thickness, porosity and hydrate saturation. After kriging interpolation of well locations, the area of each grid is almost the same, so the effective thickness, porosity and hydrate saturation are evaluated comprehensively. The background of the TOPSIS algorithm is the practical need to make a choice of advantages and disadvantages for multi-objective and multi-attribute decision-making problems, and the basic principle of the TOPSIS algorithm is to evaluate the advantages and disadvantages of the studied objectives by ranking them according to their proximity to the ideal objectives. The first step in the calculation of this method is to construct the weighted specification matrix and use the Euclidean algorithm to determine the optimal solution and the worst solution of different indicators, and then rank the evaluated objects according to the size of their distance from the optimal solution and the worst solution. If the distance of the evaluated object to the optimal solution is the smallest, and the distance to the worst solution is the largest, then the evaluated object is the best value and vice versa, and the other evaluated objects need to be ranked according to their distance from the optimal value and the size of the worst value. The calculation steps of the TOPSIS algorithm are as follows: Step 1: Establish the evaluation matrix. With n evaluation objects and m evaluation indicators, establish the original evaluation matrix X . denotes the value corresponding to the j-th evaluation object under the i- th evaluation indicator. In this paper, the evaluation matrix is the feature matrix after interpolation obtained in the previous section. X ⋯ ⋯ ⋮ ⋮ ⋮ ⋮ ⋯ (12) Step 2: Normalisation of the original matrix to eliminate the effects of differences in orders of magnitude and to harmonise the scales and units between evaluation indicators. If the evaluation indicators are positive, i.e., bigger is better, then there are: y X X X X (13) If the evaluation indicator is an inverse indicator, i.e., smaller is better, then there is: (14) If the evaluation indicator is a moderateness indicator, the optimal value is , i.e., it is best to take the value in a certain interval then there is: , , (15) At the same time, because some of the evaluation index values remain negative after the dimensionless processing, 11 leading to the meaninglessness of the logarithmic value taken when using the entropy value method to find the weights, it is necessary to non-negativity of the evaluation indexes to obtain a new matrix . When there exists 0, the translation coordinates: , 0.0001 (16) Step 3: Normalisation of the normalisation matrix. P y / y (17) Step 4: Construct the normalised decision matrix. The calculation formula is as follows: / (18) Obtain the normalised matrix : X X X ⋯ X X X ⋯ X ⋮ ⋮ ⋮ ⋮ X X ⋯ X (19) Step 5: Determine the positive and negative ideal solutions. Based on the weighted normalisation matrix, determine the positive and negative ideal sets P+ and P-. The P+ set consists of two types of indicators, for positive indicators, the maximum of the observed values is selected, and for negative indicators, the minimum of the observed values is selected; the opposite is true for the determination of the negative ideal set P-. This is expressed as follows: P X ,⋯ , X ,⋯ , X (20) P X ,⋯ , X ,⋯ , X (21) Step 6: Calculate the Euclidean distance. On the basis of determining the positive and negative ideal sets, calculate the Euclidean distance of n evaluation objects. Where Vj+ indicates the size of the distance to positive ideals and Vj- indicates the size of the distance to negative ideals, the specific formula is as follows: V X X (22) V X X (23) Step 7: Calculate relative progress Calculate the relative progress and rank the n evaluation objects to measure the amount of gas hydrate resources. / (24) Fig 4. Range of natural gas hydrate resources 12 The range of natural gas hydrate resource distribution is solved as shown in Fig. 4. As shown in Fig. 4, it can be seen that the gas hydrate resources are mainly distributed in the upper left, lower right and directly below the detection area. 5. Conclusion The results of this study demonstrate that natural gas hydrate resources are predominantly concentrated in specific areas within the exploration region. By establishing indicators and employing advanced algorithms such as Kriging and TOPSIS, we have effectively evaluated the quantity and distribution of these resources. The findings offer valuable insights for guiding future exploration and development activities related to natural gas hydrates. This research contributes to the ongoing efforts to enhance the exploration and evaluation system for natural gas hydrate resources. References [1] X.-H. Liu et al., “Evaluation of natural gas hydrate resources in the South China Sea using a new genetic analogy method,” Petroleum Science, vol. 19, no. 1, pp. 48–57, Feb. 2022. [2] X.-W. Zhang et al., “Evaluation of natural gas hydrate resources in the South China Sea by combining volumetric and trend-analysis methods,” Petroleum Science, vol. 19, no. 1, pp. 37–47, Feb. 2022. [3] Q. Meng et al., “Growth behavior and resource potential evaluation of gas hydrate in core fractures in Qilian Mountain permafrost area, Qinghai-Tibet Plateau,” China Geology, vol. 6, no. 2, pp. 208–215, Apr. 2023. [4] E. Wang, Y. Feng, T. Guo, and M. Li, “Oil content and resource quality evaluation methods for lacustrine shale: A review and a novel three-dimensional quality evaluation model,” Earth- Science Reviews, vol. 232, p. 104134, Sep. 2022. [5] J. Li et al., “Evaluating the total oil yield using a single routine Rock-Eval experiment on as-received shales,” Journal of Analytical and Applied Pyrolysis, vol. 144, p. 104707, Nov. 2019.