A novel selection method of seismic attributes based on gray relational degree and support vector machine

The selection of seismic attributes is a key process in reservoir prediction because the prediction accuracy relies on the reliability and credibility of the seismic attributes. However, effective selection method for useful seismic attributes is still a challenge. This paper presents a novel selection method of seismic attributes for reservoir prediction based on the gray relational degree (GRD) and support vector machine (SVM). The proposed method has a two-hierarchical structure. In the first hierarchy, the primary selection of seismic attributes is achieved by calculating the GRD between seismic attributes and reservoir parameters, and the GRD between the seismic attributes. The principle of the primary selection is that these seismic attributes with higher GRD to the reservoir parameters will have smaller GRD between themselves as compared to those with lower GRD to the reservoir parameters. Then the SVM is employed in the second hierarchy to perform an interactive error verification using training samples for the purpose of determining the final seismic attributes. A real-world case study was conducted to evaluate the proposed GRD-SVM method. Reliable seismic attributes were selected to predict the coalbed methane (CBM) content in southern Qinshui basin, China. In the analysis, the instantaneous amplitude, instantaneous bandwidth, instantaneous frequency, and minimum negative curvature were selected, and the predicted CBM content was fundamentally consistent with the measured CBM content. This real-world case study demonstrates that the proposed method is able to effectively select seismic attributes, and improve the prediction accuracy. Thus, the proposed GRD-SVM method can be used for the selection of seismic attributes in practice.


Introduction
Seismic attributes have been used an important feature for the purpose of reservoir prediction to recognize reservoir patterns [1][2][3][4]. However, the relationships between seismic attributes a1111111111 a1111111111 a1111111111 a1111111111 a1111111111 and reservoir lithology, fluid properties, and other parameters are rather complex. Several types of seismic attributes can even cause various adverse effects for reservoir prediction. For example, the instantaneous phase and frequency attributes are important information in seismic exploration, but these seismic attributes cannot show the stratigraphic boundary clearly [5]. The coherence attribute can increase the geological information on the time slice, but it damages the information on the vertical section at the same time [6]. Kalkomey found that the introduction of irrelevant attributes for reservoir prediction can lead to false predictions [7]. Barnes suggested that redundant seismic attributes in the candidate pool may confuse seismic interpretation [8]. Redundant seismic attributes also occupy a large storage space and require a high computational time [9]. In addition, the interaction and inner relationship between different seismic attributes in a large attributes space may lead to repetition and waste of resources. Hence, it is imperative to select the most useful/reliable seismic attributes to perform precise and efficient reservoir prediction [10][11][12].
In order for the selection of seismic attributes, Hampson et al. proposed the stepwise regression method [13]. Dorrington et al. presented the selection method of seismic attribute by using a genetic algorithm, which can select the optimal number and type of seismic attributes for the prediction of porosity [10]. Gao et al. presented a novel selection method called constrained main component analysis [14]. Ahmed et al. introduced the abductive networks to predict reservoir properties from seismic attributes. The abductive networks simultaneously selected the most relevant attributes and constructed an optimal nonlinear predictor [12]. Zhang et al. proposed a selection method of seismic attributes based on SVM, which selected attributes from various types of seismic attributes [15]. Iturrarán-Viveros used the Gamma test to analyze the data, selected the combination of seismic attributes that have the smaller Gamma statistic and trained the Artificial Neural Network (ANN) to learn to estimate porosity [16]. Qi et al. selected nine attributes that have high coefficients with the measured coalbed methane (CBM) content value, and then performed correlation analysis among these nine attributes and selected four seismic attributes with the smallest correlation coefficients to identify the CBM-enriched areas [17]. Wang et al. combined rough sets and K-L transform approach to select seismic attributes [4]. Galvis et al. presented a methodology that uses pattern recognition to select the best seismic attributes [18]. Gholamiet al. used a committee model with bat-inspired optimization algorithm to estimate porosity based on seismic attributes [19]. These excellent researches have addressed the problem of seismic attribute selection by providing the solutions from different viewpoint. However, to the best of our knowledge, the relationships between the seismic attributes and between the attributes and reservoir parameters have not been exploited in the existing selection process yet. It is worth investigating these relationships in the attribute selection because existing publications have indicated close connections between attributes and reservoir parameters [13][14][15][16][17][18][19] and can be discovered by advanced signal processing methods [20][21][22].
In order to select the most useful seismic attributes, this paper presents a novel selection method by addressing the relationships between the seismic attributes and reservoir parameters. A two-hierarchical structure was established based on the integration of gray relational degree (GRD) and support vector machine (SVM). A real world case study in southern Qinshui basin, China, has been carried out to examine the performance of the proposed GRD-SVM method. The analysis results show satisfactory prediction accuracy on the CBM content.
The rest of this paper is organized as follows. Firstly, the basic principles of GRD and SVM are elaborated, and then the basic idea of selection of seismic attributes based on GRD-SVM is introduced. Thirdly, the case study using the real world data is carried out to illustrate the effectiveness of the proposed method. Lastly, the conclusions are drawn.

Gray relational degree (GRD)
GRD was pioneered by Deng in 1984 [23]. GRD provides a simple scheme to analyze the relationships of time series, even the information provided is little informative about the objective. GRD is an analysis tool in essence, which is based on the distance between reference and comparative positions by its geometry characteristics. For a discrete time series, if the increment of two discrete time series is equal or close to equal in one period of time, the GRD of the two series is high; otherwise, the GRD is low. The basic calculation steps of GRD are as follows: Dt k Dt k \Dt kÀ 1 ¼ ;, k = 2,3,. . .,n. The two discrete time series in the time interval [a,b] are X 1 = (x 1 (t 1 ),x 1 (t 2 ),Á Á Á,x 1 (t n )), . .,n) represents the increment of the time series in time t k−1 to t k .
. . . ; nÞ represents the equalization value of increment in time t k−1 to t k of the time series. The gray correlation coefficient ξ of time series X 1 and X 2 in time t k−1 to t k is provided below: Where k is a serial number of time. Eq (1) is used to calculate the gray correlation coefficient of time series X 1 and X 2 in time t k−1 to t k . sgn(z 1 (t k ).z 2 (t k )) is a symbol function, which reflects the correlations of the two sequences. In addition, when z 1 (t k ).z 2 (t k )>0, then ξ(t k )>0, and when z 1 (t k ).z 2 (t k )<0, then ξ(t k )<0.kz 1 (t k )|−|z 2 (t k )k is the increment of the constitute difference, and ð1 À minðjz 1 ðt k Þj;jz 2 ðt k ÞjÞ maxðjz 1 ðt k Þj;jz 2 ðt k ÞjÞ Þ is the constituent ratio. When the increment between the two time series in time t k-1 to t k is equal or close to equal, kz 1 (t k )|−|z 2 (t k )k and ð1 À minðjz 1 ðt k Þj;jz 2 ðt k ÞjÞ maxðjz 1 ðt k Þj;jz 2 ðt k ÞjÞ Þ are close to zero. At this time, the gray correlation coefficient of the two time series in time t k−1 to t k is high (close to 1); otherwise, the gray correlation coefficient is low.
After we calculate the gray correlation coefficient, the GRD in time interval [a,b] can be calculated by the Eq (2).
Where r is the GRD of the time series X 1 and X 2 . When we assume t k = k, k = 1,2,. . .,n, Eq (2) can be transformed to the following constrained formation:

Support vector machine (SVM)
SVM, which is a set of related supervised learning methods for classification and regression, was first introduced by Vapnik [24]. SVM provides the best generalization accuracy classifier by maximizing the margin between two classes in a feature space. Because of this, it has been widely applied in several fields such as classification, regression, face recognition, and time series prediction. Support vector regression is an application of SVM for function regression. The basic theory of SVM has been detailed described in several literatures, and thus, the description of the theory will not be repeated in this paper. We consider linear support vector regression machine as an example and briefly introduce the basic principle of SVM.
In a given training sample set {(x i ,y i ),i = 1,2. . .,l}, x i 2 R n represents the input value, y i 2 R represents the corresponding target value, l is the number of samples. f(x) = <w,x>+b is used for fitting the sample data set {x i ,y i }, and b is the offset value.
In Eq (4), ε is the insensitive loss function, which is directly related to function estimation precision, and its value is greater than zero. The main purpose of this study is to construct a function f(x), and, at the same time, ensure that the distance between the function and target is less than ε. Thus, for those unknown samples x, the function f(x) can optimally estimate the corresponding target value.
In the linear case, we can assume f(x) for the following form: Where w 2 R n is the weight vector of function f(x). x 2 R n is the input value, <.> is the inner product, and b 2 R is the offset value of function f(x).
For a given training sample, (x i ,y i ), i = 1,2,. . .,l, x 2 R n , y 2 R, the optimization problem is as follows: Subject to: l When a given constraint condition is not able to solve Eq (6), the slack variables ξ i and x Ã i are introduced, which represent the distance from the actual values to the corresponding boundary values of ε (the insensitive loss function). Eq (6) can be transformed into the following constrained formation: Subject to: l Where C is the penalty coefficient, which is taken as the regularized constant that determines the trade-off between the empirical error (risk) and regularization term. Eq (7) represents quadratic programming problems with linear inequality constraints, which can be solved using the Lagrange multiplier method, i.e., Where a i and a Ã i are called Lagrange multipliers.

Selection method of seismic attributes based on GRD and SVM
GRD can truly reflect the degree of closeness of relative changes in the time series. The GRD between seismic attributes and reservoir parameters, and the GRD between different seismic attributes can be calculated to achieve primary selection of seismic attributes. SVM is very suitable for small samples, nonlinear cases, and other issues. It has better generalization and promotion capabilities. Thus, SVM can be used to perform interactive error verification of known samples and achieve final selection of seismic attributes. The flow chart of the selection method of seismic attributes based on GRD and SVM is shown in Fig 1. The basic idea of GRD-SVM is as follows: 1. GRD is used to select seismic attributes for preprocessing, i.e., based on the GRD between seismic attributes and reservoir parameters, 6-8 types of seismic attributes are initially selected, which have high correlation with reservoir parameters, in order to effectively remove redundant attributes and reduce the number of seismic attributes used for reservoir prediction.
2. The 6-8 types of seismic attributes that were selected may have high GRD with each other. The GRD between the selected seismic attributes is calculated, and seismic attributes with higher GRDs are clustered into one category. The seismic attribute that has the largest GRD with reservoir parameters is selected from each category; thus, 3-5 types of sensitive seismic attributes are selected.
3. SVM is used to perform interactive error verification for known samples. When a seismic attribute is selected, a particular sample is excluded from testing. The rest of the samples in the training set are used to calculate the predicted root mean square error (RMSE). The above procedure is repeated for all the samples. The calculated RMSEs of all the samples are added and then divided by the total number of samples to obtain the RMSE of one seismic attribute. The number of seismic attributes is increased, and the above steps are repeated. The RMSE of different types of seismic attributes are calculated, and seismic attributes with the smallest RMSE are selected.

Case study
CBM is one type of gas spontaneously occurring in coal seams. The prediction of CBM content is a crucial factor for production safety in coal mines and formulation of plans for the development of CBM resources. The use of seismic attributes to predict CBM content has made great progress [25,26], but there is a lack of research on selection methods of seismic attributes for the prediction of CBM content. The research area was located at the southern part of the Qinshui Basin. The coal series was developed mainly in the Lower Permian Shanxi group and the Upper Carboniferous Taiyuan group. A total of 15 coal beds were found in these strata with a mean total thickness 136.02 m. Three coal seams in the Shanxi Group and fifteen coal seams in the Taiyuan Group are minable. The thickness of the three coal seams in the Shanxi Group ranges from 6.49 to 7.45 m with an average thickness of 6.79 m. The main structure is typical of synclinal and anticlinal composite folds extending to NNE and NE, respectively [17].
Qi and Zhang [17] studied the seismic attributes of this area. A total of 64 seismic attributes, including coal seam reflected waves, amplitude, complex seismic trace, and sequence statistics attributes were used. Table 1 is the seismic attributes used in this study, the calculation method of each seismic attribute can see the help manual of Landmark and Geomodeling. The selected seismic attributes were combined based on the Dempster-Shafer evidence theory, and CBM enriched areas were predicted. In order to verify the selection method of seismic attributes proposed in this paper, the seismic attributes of Qi's literature were used to select and analyze.
The measured CBM content of the seven wells located in 3# coal seam of Qinshui Basin were 18.02m 3 /t, 17.58m 3 /t, 16.86m 3 /t, 12.51m 3 /t, 10.12m 3 /t, 9.79m 3 /t and 8.68m 3 /t, respectively. Based on the value of GRD between the seismic attributes and CBM content, we selected instantaneous amplitude, instantaneous bandwidth, instantaneous main frequency, instantaneous frequency, minimum negative curvature, instantaneous phase, maximum positive curvature, and absorption attribute. Then we calculated the GRDs between the eight types of seismic attributes, which were shown in Table 2.
As shown in Table 2, some seismic attributes have high GRD with other seismic attributes. The GRD between instantaneous bandwidth and maximum positive curvature is 0.7049, and the GRD between instantaneous amplitude and instantaneous main frequency is -0.5797. Based on the selection principle of seismic attributes with a greater GRD, we further selected instantaneous amplitude, instantaneous bandwidth, instantaneous frequency, minimum negative curvature, instantaneous phase, and absorption attribute as sensitive seismic attributes. SVM was used to perform interactive error verification between the sensitive seismic attributes and the known CBM content, and the final selection of seismic attributes was achieved based on the RMSE between the predicted and actual values. The specific process is as follows: 1. When a particular type of seismic attribute was selected, the data of one well was removed from the training data set (this well is usually known as the blind well). The SVM training data of the other six wells were used, and the obtained training results were used to predict the result of the blind well. The RMSE between the predicted and the theoretical values was calculated.
2. This process was sequentially performed for all the seven wells present in this work area. The RMSE of the seismic attribute was added and then divided by 7 to obtain the RMSE of this seismic attribute.
3. The types of seismic attributes were increased, and the above process was repeated to obtain the RMSE of different types of seismic attributes, as shown in Table 3. As shown in Table 3, the minimum RMSE between the predicted and actual values is 0.6305, when four seismic attributes are selected. Thus, the final seismic attributes selected are as follows: instantaneous amplitude, instantaneous frequency, instantaneous bandwidth, and minimum negative curvature. The four seismic attributes selected along 3coal seam of Qinshui Basin is shown in Figs 2-5. Wells with a high CBM content are located in the strong instantaneous amplitude area, and wells with a low CBM content are located in the weak instantaneous amplitude area.
(The marked values in the figure is content of CBM for wells, the unit is m 3 /t) As shown in Figs 3-5, the well with a high CBM content is located in the area with a low instantaneous bandwidth, low instantaneous frequency, and high minimum negative curvature, whereas the well with a low CBM content is located in the area with a high instantaneous bandwidth, high instantaneous frequency, and low minimum negative curvature. However, the correlation between the CBM content and these three seismic attributes is not as significant as instantaneous amplitude.
The four types of seismic attributes were used to predict the CBM content of the whole work area based on SVM, and the results are shown in Fig 6.  Fig 6 shows the prediction results of CBM content. In Fig 6, two wells (measured CBM con-tent18.02m 3 /t and 16.85m 3 /t) are located in areas with a high CBM content. Three wells (measured CBM content 10.12m 3 /t, 9.79m 3 /t, and 8.68m 3 /t) are located in areas with the low CBM content. Two wells (measured CBM content 17.58m 3 /t and 12.51m 3 /t) are located in areas with medium CBM content. The prediction results were compared with the measured CBM content, and there is only one well (measured CBM content 17.58m 3 /t) that does not match with the actual value. The prediction results are mainly consistent with the actual values with an accuracy of 85.7%.
A certain thickness of coal seam is the base of CBM reservoir. The coal reservoir provides both sources and reservoir spaces for CBM. The thicker the coal seam, the better the production of CBM [27].   The area possible of containing high CBM is between the well with a measured CBM content of 12.51m 3 /t and the well with a measured CBM content of 16.85m 3 /t. The area possible of containing high CBM (Fig 9) is different from the area of the thick coal seam (Fig 8). However, there are three areas with high CBM contents near the location of the wells (measured CBM contents of 16.85m 3 /t, 12.51m 3 /t, and 18.02m 3 /t), as shown in Fig 6. The scope and positions of areas with high CBM content, which were predicted by this study (Fig 6) are essentially consistent with the thickness of the coal seam.
Thus by contrasting and analyzing Fig 6, Fig 8 and Fig 9, the selection method of seismic attributes based on GRD and SVM could effectively select seismic attributes, and the selection results of seismic attributes are more reasonable. The scope and positions of areas with high CBM content, which were predicted by the selection of seismic attributes, are basically consistent with the thickness of the coal seam.
We also use the 8 seismic attributes: instantaneous amplitude, instantaneous bandwidth, instantaneous main frequency, instantaneous frequency, minimum negative curvature, instantaneous phase, maximum positive curvature, and absorption attribute to carry out Step-wise regression. Fig 10 is a validation plot for the Step-wise regression analysis. The blue (lower) curve shows the error calculated using the Training Data. The red (upper) curve shows the error calculated using the Validation Data. Fig 10 shows that when 5 seismic attributes are selected, the validation error is the smallest. When 6 or more attributes are used, the validation error increases, meaning that these additional attributes are over-fitting the data. The RMSE of Step-wise regression is much larger than the method proposed in this paper. This also illustrates the effectiveness of the method presented in this paper.

Conclusions
S8 File. The data is used in "application to seismic data" section. (MAT) S9 File. The data is used in "application to seismic data" section. (TXT)