notebook.community

Edit and run

Facies classification using ExtraTrees, OneVsRest and Feature imputation

Contest entry by :geoLEARN

Martin Blouin, Lorenzo Perozzi, Antoine Caté

in collaboration with Erwan Gloaguen

Original contest notebook by Brendon Hall, Enthought

In this notebook we will train a machine learning algorithm to predict facies from well log data. The dataset comes from a class exercise from The University of Kansas on Neural Networks and Fuzzy Systems. This exercise is based on a consortium project to use machine learning techniques to create a reservoir model of the largest gas fields in North America, the Hugoton and Panoma Fields. For more info on the origin of the data, see Bohling and Dubois (2003) and Dubois et al. (2007).

The dataset consists of log data from nine wells that have been labeled with a facies type based on observation of core. We will use this log data to train a Random Forest model to classify facies types.

Exploring the dataset

First, we import and examine the dataset used to train the classifier.





In [1]:



    


###### Importing all used packages
%matplotlib inline

import warnings
warnings.filterwarnings('ignore')

import pandas as pd
import numpy as np

from pandas import set_option
pd.options.mode.chained_assignment = None

###### Import packages needed for the make_vars functions
import Feature_Engineering as FE
import xgboost as xgb
##### import stuff from scikit learn
from sklearn.ensemble import ExtraTreesClassifier, RandomForestRegressor
from sklearn.multiclass import OneVsRestClassifier
from sklearn.model_selection import LeaveOneGroupOut, cross_val_predict

filename = 'facies_vectors.csv'
df = pd.read_csv(filename)
df.head()



    


















    






/usr/local/anaconda3/lib/python3.5/site-packages/sklearn/cross_validation.py:44: DeprecationWarning: This module was deprecated in version 0.18 in favor of the model_selection module into which all the refactored classes and functions are moved. Also note that the interface of the new CV iterators are different from that of this module. This module will be removed in 0.20.
  "This module will be removed in 0.20.", DeprecationWarning)










    
Out[1]:










  
    

      
      Facies
      Formation
      Well Name
      Depth
      GR
      ILD_log10
      DeltaPHI
      PHIND
      PE
      NM_M
      RELPOS
    

  
  
    

      0
      3
      A1 SH
      SHRIMPLIN
      2793.0
      77.45
      0.664
      9.9
      11.915
      4.6
      1
      1.000
    

    

      1
      3
      A1 SH
      SHRIMPLIN
      2793.5
      78.26
      0.661
      14.2
      12.565
      4.1
      1
      0.979
    

    

      2
      3
      A1 SH
      SHRIMPLIN
      2794.0
      79.05
      0.658
      14.8
      13.050
      3.6
      1
      0.957
    

    

      3
      3
      A1 SH
      SHRIMPLIN
      2794.5
      86.10
      0.655
      13.9
      13.115
      3.5
      1
      0.936
    

    

      4
      3
      A1 SH
      SHRIMPLIN
      2795.0
      74.58
      0.647
      13.5
      13.300
      3.4
      1
      0.915
    

  




















In [2]:



    


####### create X_train and y_train
X_train = df[['GR', 'ILD_log10', 'DeltaPHI', 'PHIND', 'NM_M', 'RELPOS']][df.PE.notnull()]
y_train = df['PE'][df.PE.notnull()]
groups_train = df['Well Name'][df.PE.notnull()]

X_fit = df[['GR', 'ILD_log10', 'DeltaPHI', 'PHIND', 'NM_M', 'RELPOS']][df.PE.isnull()]



    











In [3]:



    


Cl = RandomForestRegressor(n_estimators=100)
Cl.fit(X_train, y_train)
y_predict = Cl.predict(X_fit)
df['PE'][df.PE.isnull()] = y_predict
training_data = df

Feature engineering

As stated in our previous submission, we believe that feature engineering has a high potential for increasing classification success. A strategy for building new variables is explained below.

The dataset is distributed along a series of drillholes intersecting a stratigraphic sequence. Sedimentary facies tend to be deposited in sequences that reflect the evolution of the paleo-environment (variations in water depth, water temperature, biological activity, currents strenght, detrital input, ...). Each facies represents a specific depositional environment and is in contact with facies that represent a progressive transition to an other environment. Thus, there is a relationship between neighbouring samples, and the distribution of the data along drillholes can be as important as data values for predicting facies.

A series of new variables (features) are calculated and tested below to help represent the relationship of neighbouring samples and the overall texture of the data along drillholes. These variables are:

  • detail and approximation coeficients at various levels of two wavelet transforms (using two types of Daubechies wavelets);
  • measures of the local entropy with variable observation windows;
  • measures of the local gradient with variable observation windows;
  • rolling statistical calculations (i.e., mean, standard deviation, min and max) with variable observation windows;
  • ratios between marine and non-marine lithofacies with different observation windows;
  • distances from the nearest marine or non-marine occurence uphole and downhole.

Functions used to build these variables are located in the Feature Engineering python script.

All the data exploration work related to the conception and study of these variables is not presented here.





In [4]:



    


##### cD From wavelet db1
dwt_db1_cD_df = FE.make_dwt_vars_cD(wells_df=training_data, logs=['GR', 'ILD_log10', 'DeltaPHI', 'PE', 'PHIND'],
                       levels=[1, 2, 3, 4], wavelet='db1')

##### cA From wavelet db1
dwt_db1_cA_df = FE.make_dwt_vars_cA(wells_df=training_data, logs=['GR', 'ILD_log10', 'DeltaPHI', 'PE', 'PHIND'],
                       levels=[1, 2, 3, 4], wavelet='db1')

##### cD From wavelet db3
dwt_db3_cD_df = FE.make_dwt_vars_cD(wells_df=training_data, logs=['GR', 'ILD_log10', 'DeltaPHI', 'PE', 'PHIND'],
                       levels=[1, 2, 3, 4], wavelet='db3')

##### cA From wavelet db3
dwt_db3_cA_df = FE.make_dwt_vars_cA(wells_df=training_data, logs=['GR', 'ILD_log10', 'DeltaPHI', 'PE', 'PHIND'],
                       levels=[1, 2, 3, 4], wavelet='db3')

##### From entropy
entropy_df = FE.make_entropy_vars(wells_df=training_data, logs=['GR', 'ILD_log10', 'DeltaPHI', 'PE', 'PHIND'],
                               l_foots=[2, 3, 4, 5, 7, 10])

###### From gradient
gradient_df = FE.make_gradient_vars(wells_df=training_data, logs=['GR', 'ILD_log10', 'DeltaPHI', 'PE', 'PHIND'],
                                 dx_list=[2, 3, 4, 5, 6, 10, 20])

##### From rolling average
moving_av_df = FE.make_moving_av_vars(wells_df=training_data, logs=['GR', 'ILD_log10', 'DeltaPHI', 'PE', 'PHIND'],
                                   windows=[1, 2, 5, 10, 20])

##### From rolling standard deviation
moving_std_df = FE.make_moving_std_vars(wells_df=training_data, logs=['GR', 'ILD_log10', 'DeltaPHI', 'PE', 'PHIND'],
                                     windows=[3 , 4, 5, 7, 10, 15, 20])

##### From rolling max
moving_max_df = FE.make_moving_max_vars(wells_df=training_data, logs=['GR', 'ILD_log10', 'DeltaPHI', 'PE', 'PHIND'],
                                     windows=[3, 4, 5, 7, 10, 15, 20])

##### From rolling min
moving_min_df = FE.make_moving_min_vars(wells_df=training_data, logs=['GR', 'ILD_log10', 'DeltaPHI', 'PE', 'PHIND'],
                                     windows=[3 , 4, 5, 7, 10, 15, 20])

###### From rolling NM/M ratio
rolling_marine_ratio_df = FE.make_rolling_marine_ratio_vars(wells_df=training_data, windows=[5, 10, 15, 20, 30, 50, 75, 100, 200])

###### From distance to NM and M, up and down
dist_M_up_df = FE.make_distance_to_M_up_vars(wells_df=training_data)
dist_M_down_df = FE.make_distance_to_M_down_vars(wells_df=training_data)
dist_NM_up_df = FE.make_distance_to_NM_up_vars(wells_df=training_data)
dist_NM_down_df = FE.make_distance_to_NM_down_vars(wells_df=training_data)



    











In [5]:



    


list_df_var = [dwt_db1_cD_df, dwt_db1_cA_df, dwt_db3_cD_df, dwt_db3_cA_df,
               entropy_df, gradient_df, moving_av_df, moving_std_df, moving_max_df, moving_min_df,
              rolling_marine_ratio_df, dist_M_up_df, dist_M_down_df, dist_NM_up_df, dist_NM_down_df]
combined_df = training_data
for var_df in list_df_var:
    temp_df = var_df
    combined_df = pd.concat([combined_df,temp_df],axis=1)
combined_df.replace(to_replace=np.nan, value='-1', inplace=True)
print (combined_df.shape)
combined_df.head(5)



    


















    






(4149, 299)










    
Out[5]:










  
    

      
      Facies
      Formation
      Well Name
      Depth
      GR
      ILD_log10
      DeltaPHI
      PHIND
      PE
      NM_M
      ...
      Marine_ratio_20_centered
      Marine_ratio_30_centered
      Marine_ratio_50_centered
      Marine_ratio_75_centered
      Marine_ratio_100_centered
      Marine_ratio_200_centered
      dist_M_up
      dist_M_down
      dist_NM_up
      dist_NM_down
    

  
  
    

      0
      3
      A1 SH
      SHRIMPLIN
      2793.0
      77.45
      0.664
      9.9
      11.915
      4.6
      1
      ...
      1.0
      1.0
      1.0
      1.0
      1.140000
      1.510000
      -1.0
      21.5
      0.0
      0.0
    

    

      1
      3
      A1 SH
      SHRIMPLIN
      2793.5
      78.26
      0.661
      14.2
      12.565
      4.1
      1
      ...
      1.0
      1.0
      1.0
      1.0
      1.156863
      1.504950
      -1.0
      21.0
      0.0
      0.0
    

    

      2
      3
      A1 SH
      SHRIMPLIN
      2794.0
      79.05
      0.658
      14.8
      13.050
      3.6
      1
      ...
      1.0
      1.0
      1.0
      1.0
      1.173077
      1.500000
      -1.0
      20.5
      0.0
      0.0
    

    

      3
      3
      A1 SH
      SHRIMPLIN
      2794.5
      86.10
      0.655
      13.9
      13.115
      3.5
      1
      ...
      1.0
      1.0
      1.0
      1.0
      1.188679
      1.495146
      -1.0
      20.0
      0.0
      0.0
    

    

      4
      3
      A1 SH
      SHRIMPLIN
      2795.0
      74.58
      0.647
      13.5
      13.300
      3.4
      1
      ...
      1.0
      1.0
      1.0
      1.0
      1.203704
      1.490385
      -1.0
      19.5
      0.0
      0.0
    

  



5 rows × 299 columns

Building a prediction model from these variables





In [21]:



    


X = combined_df.iloc[:, 4:]
y = combined_df['Facies']
groups = combined_df['Well Name']
XGBoost = xgb.XGBClassifier(n_estimators=200, max_depth=5, min_child_weight = 2, learning_rate=0.12, subsample=0.7, colsample_bytree=0.9, seed=i)
Cl = OneVsRestClassifier(XGBoost, n_jobs=-1)

Applying the classification model to test data





In [22]:



    


filename = 'validation_data_nofacies.csv'
test_data = pd.read_csv(filename)
test_data.head(5)



    


















    
Out[22]:










  
    

      
      Formation
      Well Name
      Depth
      GR
      ILD_log10
      DeltaPHI
      PHIND
      PE
      NM_M
      RELPOS
    

  
  
    

      0
      A1 SH
      STUART
      2808.0
      66.276
      0.630
      3.3
      10.65
      3.591
      1
      1.000
    

    

      1
      A1 SH
      STUART
      2808.5
      77.252
      0.585
      6.5
      11.95
      3.341
      1
      0.978
    

    

      2
      A1 SH
      STUART
      2809.0
      82.899
      0.566
      9.4
      13.60
      3.064
      1
      0.956
    

    

      3
      A1 SH
      STUART
      2809.5
      80.671
      0.593
      9.5
      13.25
      2.977
      1
      0.933
    

    

      4
      A1 SH
      STUART
      2810.0
      75.971
      0.638
      8.7
      12.35
      3.020
      1
      0.911
    

  




















In [23]:



    


##### cD From wavelet db1
dwt_db1_cD_df = FE.make_dwt_vars_cD(wells_df=test_data, logs=['GR', 'ILD_log10', 'DeltaPHI', 'PE', 'PHIND'],
                                    levels=[1, 2, 3, 4], wavelet='db1')

##### cA From wavelet db1
dwt_db1_cA_df = FE.make_dwt_vars_cA(wells_df=test_data, logs=['GR', 'ILD_log10', 'DeltaPHI', 'PE', 'PHIND'],
                       levels=[1, 2, 3, 4], wavelet='db1')

##### cD From wavelet db3
dwt_db3_cD_df = FE.make_dwt_vars_cD(wells_df=test_data, logs=['GR', 'ILD_log10', 'DeltaPHI', 'PE', 'PHIND'],
                       levels=[1, 2, 3, 4], wavelet='db3')

##### cA From wavelet db3
dwt_db3_cA_df = FE.make_dwt_vars_cA(wells_df=test_data, logs=['GR', 'ILD_log10', 'DeltaPHI', 'PE', 'PHIND'],
                       levels=[1, 2, 3, 4], wavelet='db3')

##### From entropy
entropy_df = FE.make_entropy_vars(wells_df=test_data, logs=['GR', 'ILD_log10', 'DeltaPHI', 'PE', 'PHIND'],
                               l_foots=[2, 3, 4, 5, 7, 10])

###### From gradient
gradient_df = FE.make_gradient_vars(wells_df=test_data, logs=['GR', 'ILD_log10', 'DeltaPHI', 'PE', 'PHIND'],
                                 dx_list=[2, 3, 4, 5, 6, 10, 20])

##### From rolling average
moving_av_df = FE.make_moving_av_vars(wells_df=test_data, logs=['GR', 'ILD_log10', 'DeltaPHI', 'PE', 'PHIND'],
                                   windows=[1, 2, 5, 10, 20])

##### From rolling standard deviation
moving_std_df = FE.make_moving_std_vars(wells_df=test_data, logs=['GR', 'ILD_log10', 'DeltaPHI', 'PE', 'PHIND'],
                                     windows=[3 , 4, 5, 7, 10, 15, 20])

##### From rolling max
moving_max_df = FE.make_moving_max_vars(wells_df=test_data, logs=['GR', 'ILD_log10', 'DeltaPHI', 'PE', 'PHIND'],
                                     windows=[3, 4, 5, 7, 10, 15, 20])

##### From rolling min
moving_min_df = FE.make_moving_min_vars(wells_df=test_data, logs=['GR', 'ILD_log10', 'DeltaPHI', 'PE', 'PHIND'],
                                     windows=[3 , 4, 5, 7, 10, 15, 20])

###### From rolling NM/M ratio
rolling_marine_ratio_df = FE.make_rolling_marine_ratio_vars(wells_df=test_data, windows=[5, 10, 15, 20, 30, 50, 75, 100, 200])

###### From distance to NM and M, up and down
dist_M_up_df = FE.make_distance_to_M_up_vars(wells_df=test_data)
dist_M_down_df = FE.make_distance_to_M_down_vars(wells_df=test_data)
dist_NM_up_df = FE.make_distance_to_NM_up_vars(wells_df=test_data)
dist_NM_down_df = FE.make_distance_to_NM_down_vars(wells_df=test_data)



    











In [24]:



    


combined_test_df = test_data
list_df_var = [dwt_db1_cD_df, dwt_db1_cA_df, dwt_db3_cD_df, dwt_db3_cA_df,
               entropy_df, gradient_df, moving_av_df, moving_std_df, moving_max_df, moving_min_df,
              rolling_marine_ratio_df, dist_M_up_df, dist_M_down_df, dist_NM_up_df, dist_NM_down_df]
for var_df in list_df_var:
    temp_df = var_df
    combined_test_df = pd.concat([combined_test_df,temp_df],axis=1)
combined_test_df.replace(to_replace=np.nan, value='-99999', inplace=True)

X_test = combined_test_df.iloc[:, 3:]



    











In [29]:



    


test_pred_df = combined_test_df[['Well Name', 'Depth']]
for i in range(20):
    XGBoost = xgb.XGBClassifier(n_estimators=200, max_depth=5, min_child_weight = 2, learning_rate=0.12, subsample=0.7, colsample_bytree=0.9, seed=i)
    Cl = OneVsRestClassifier(XGBoost, n_jobs=-1)
    Cl.fit(X.values, y.values)
    y_test = Cl.predict(X_test.values)
    y_test = pd.DataFrame(y_test, columns=['Predicted Facies #' + str(i)])
    test_pred_df = pd.concat([test_pred_df, y_test], axis=1)
print (test_pred_df.shape)
test_pred_df.head()



    


















    






(830, 22)










    
Out[29]:










  
    

      
      Well Name
      Depth
      Predicted Facies #0
      Predicted Facies #1
      Predicted Facies #2
      Predicted Facies #3
      Predicted Facies #4
      Predicted Facies #5
      Predicted Facies #6
      Predicted Facies #7
      ...
      Predicted Facies #10
      Predicted Facies #11
      Predicted Facies #12
      Predicted Facies #13
      Predicted Facies #14
      Predicted Facies #15
      Predicted Facies #16
      Predicted Facies #17
      Predicted Facies #18
      Predicted Facies #19
    

  
  
    

      0
      STUART
      2808.0
      3
      2
      3
      3
      2
      3
      2
      2
      ...
      3
      3
      2
      3
      2
      2
      3
      2
      3
      3
    

    

      1
      STUART
      2808.5
      3
      3
      3
      3
      3
      3
      3
      3
      ...
      3
      3
      3
      3
      3
      3
      3
      3
      3
      3
    

    

      2
      STUART
      2809.0
      3
      3
      3
      3
      3
      3
      3
      3
      ...
      3
      3
      3
      3
      3
      3
      3
      3
      3
      3
    

    

      3
      STUART
      2809.5
      3
      3
      3
      3
      3
      3
      3
      3
      ...
      3
      3
      3
      3
      3
      3
      3
      3
      3
      3
    

    

      4
      STUART
      2810.0
      3
      3
      3
      3
      3
      3
      3
      3
      ...
      3
      3
      3
      3
      3
      3
      3
      3
      3
      3
    

  



5 rows × 22 columns

Exporting results





In [ ]:



    


test_pred_df.to_pickle('Prediction_submission4_XGBoost2.pkl')

Content source: seg/2016-ml-contest

Similar notebooks:

notebook.community | gallery | about