Utilizzando pykrige e matplotlib per la visualizzazione spaziale delle variazioni geologiche

Using pykrige and matplotlib for spatial visualization of geological variations.

Esplorare la Variabilità Geologica Spaziale dai Dati dei Log di Pozzo

Variabilità spaziale delle misurazioni di compressione acustica della lente sottomarina norvegese. Immagine dell'autore.

Quando si lavora con dati geologici e petrofisici, spesso vogliamo capire come quei dati cambiano nel nostro campo o nella nostra regione. Uno dei modi in cui possiamo farlo è quello di grigliare i nostri valori di misurazione effettivi e di estrapolare quali potrebbero essere quei valori in altre aree che devono ancora essere esplorate utilizzando fori di sondaggio.

Un particolare metodo per effettuare questa estrapolazione è il kriging, una procedura geostatistica denominata in onore dell’ingegnere minerario sudafricano Danie G. Krige. L’idea alla base del kriging risiede nella sua tecnica di stima: utilizza la correlazione spaziale tra i dati osservati per prevedere i valori in località non misurate.

Misurando come le variabili cambiano nel tempo, questo metodo stabilisce una relazione statistica che può essere utilizzata per prevedere i valori in un’area, trasformando punti dati sparsi in una mappa spaziale coerente.

In questo tutorial, esamineremo una libreria Python chiamata pykrige . Questa libreria è stata progettata per calcoli di kriging 2D e 3D ed è facile da usare con i dati dei log di pozzo.

Importazione di Librerie e Dati

Prima di tutto, dobbiamo importare le librerie di cui avremo bisogno. Per questo articolo, avremo bisogno delle seguenti librerie:

  • pandas — per leggere i nostri dati, che sono in formato csv
  • matplotlib per creare la nostra visualizzazione
  • pykrige per effettuare il kriging
  • numpy per alcuni calcoli numerici
import pandas as pdimport matplotlib.pyplot as pltfrom pykrige import OrdinaryKrigingimport numpy as np

Dopo aver importato le librerie, possiamo ora importare i nostri dati.

In questo tutorial, utilizzeremo un dataset derivato dalla competizione di Machine Learning Xeek e Force 2020 per la previsione della litologia dai dati dei log di pozzo. I dettagli di questo dataset possono essere trovati alla fine di questo articolo.

Questo sottoinsieme del dataset della competizione contiene 65 località di pozzo con misurazioni medie di compressione acustica per la formazione di Balder.

Per leggere i nostri dati possiamo utilizzare la funzione read_csv() di pandas e passare la posizione del file di dati. In questo esempio, utilizziamo un percorso relativo al nostro Jupyter Notebook, ma potremmo utilizzare un percorso assoluto se il nostro file si trova altrove.

df = pd.read_csv('Data/Xeek Force 2020/Xeek_2020_Balder_DTC_AVG.csv')df

Quando visualizziamo il dataframe, vediamo che abbiamo 65 pozzi, che contengono la posizione della parte superiore della formazione di Balder (coordinate di griglia X_LOC e Y_LOC e LAT & LON per la latitudine e la longitudine). Abbiamo anche la Profondità Verticale Vera Submarina (TVDSS) in cui la formazione è stata incontrata e il valore medio per la compressione acustica (DTC).

Dataframe che contiene i dati di posizione per i pozzi selezionati e il valore medio DTC - compressione acustica per la formazione di Balder nel Mare del Nord norvegese. Immagine dell'autore.

Visualizzazione delle Posizioni Spaziali dei Pozzi

Ora che i nostri dati sono stati caricati con successo in un dataframe, possiamo visualizzare i nostri dati per capire dove si trovano i nostri pozzi. Per farlo, utilizzeremo il grafico a dispersione di matplotlib e passeremo le colonne di longitudine e latitudine.

plt.scatter(df['Longitude'], df['Latitude'], c=df['DTC'])

Quando eseguiamo il codice sopra, otteniamo il seguente grafico.

Una figura di base di matplotlib mostra la posizione dei nostri pozzi e i valori DTC nella regione del Mare del Nord norvegese. Immagine dell'autore.

Possiamo vedere che la figura sopra è molto semplice, senza colorbar o etichette degli assi.

Modifichiamo leggermente il grafico aggiungendo queste funzionalità.

cm = plt.cm.get_cmap('viridis')plt.figure(figsize=(10,10))scatter = plt.scatter(df['LON'], df['LAT'], c=df['DTC_MEAN'], cmap=cm, s=50)plt.colorbar(scatter)plt.xlabel('Longitudine')plt.ylabel('Latitudine')plt.show()

Quando eseguiamo il codice sopra, otteniamo la seguente figura, che ci fornisce maggiori informazioni sui nostri dati. Possiamo utilizzare la colorbar per stimare i valori dei nostri punti.

Grafico a dispersione di Matplotlib dei pozzi selezionati nel Mare del Nord norvegese dopo l'aggiunta di colorbar e etichette. Immagine dell'autore.

Applicazione del Kriging

Per comprendere meglio i nostri punti dati e come la misura DTC varia nell’area per la formazione di Balder, possiamo utilizzare il kriging e i nostri punti dati per riempire le lacune tra i nostri valori misurati.

Per fare ciò, dobbiamo creare un oggetto OrdinaryKriging dalla libreria pykrige.

In questo oggetto passiamo i nostri dati di posizione per x e y e i dati che vogliamo mappare al parametro z.

Dobbiamo anche selezionare quale modello variogramma vogliamo utilizzare. In questo caso, utilizzeremo un modello esponenziale. Maggiori dettagli sui tipi di modelli possono essere trovati nella documentazione.

Poiché stiamo utilizzando la latitudine e la longitudine per le nostre coordinate x e y, possiamo cambiare il parametro coordinates_type in geografico

OK = OrdinaryKriging(x=df['LON'],                       y=df['LAT'],                       z=df['DTC_MEAN'],                      variogram_model='exponential',                      verbose=True, enable_plotting=True,                      coordinates_type='geographic')

Quando eseguiamo il codice sopra, otteniamo il seguente riassunto del modello e semivariogramma.

Riassunto del modello di pykrige. Immagine dell'autore.

Ecco un breve riassunto dei parametri restituiti:

  • Nugget : Il nugget è l’intercetta del variogramma, rappresentante la varianza a distanza zero, spesso dovuta a errori di misura o variazioni su piccola scala.
  • Full Sill : Il sill è la massima varianza che il variogramma raggiunge e inizia ad appiattirsi, il che accade quando i punti sono molto distanti.
  • Range : La portata è la distanza alla quale il variogramma raggiunge il sill, il che significa che la distanza oltre la quale una ulteriore separazione dei punti non aumenta la varianza.
  • Partial Sill : Il parziale sill è la differenza tra il sill e il nugget, rappresentante la quantità di varianza che è strutturata spazialmente nei dati.

Ciò può darci un’idea di quanto adatto sia il nostro modello per i dati in base alla forma della linea generata e dei punti.

Visualizzazione dei risultati del Kriging

Per iniziare a visualizzare i nostri dati, dobbiamo creare una griglia di dati.

Per farlo, creiamo prima degli array per le latitudini e le longitudini tra le coordinate che definiamo. In questo caso, vorremmo che il grafico si estendesse da 57,5 gradi N a 62 gradi N e da 1,5 gradi E a 4,5 gradi E.

Utilizzando np.arange ci permetterà di creare questi array a distanza regolare.

grid_lat = np.arange(57.5, 62, 0.01, dtype='float64')grid_long = np.arange(1.5, 4.5, 0.01,dtype='float64')

Ora che abbiamo le coordinate X & Y, possiamo creare la nostra griglia di valori. Per farlo, richiamiamo OK.execute , e passiamo i nostri array di latitudine e longitudine.

zstar, ss = OK.execute('grid', grid_long, grid_lat)

Ciò restituirà due array. La nostra griglia di dati (zstar) e l’incertezza associata ad essa (ss).

Successivamente, possiamo ora utilizzare il nostro array di dati e tracciarlo utilizzando imshow di matplotlib.

import matplotlib.pyplot as pltfig, ax = plt.subplots(figsize=(10,10))image = ax.imshow(zstar, extent=(1.5, 4.5, 57.5, 62), origin='lower')ax.set_xlabel('Longitudine', fontsize=14, fontweight='bold')ax.set_ylabel('Latitudine', fontsize=14, fontweight='bold')scatter = ax.scatter(x=df['LON'], y=df['LAT'], color='black')colorbar = fig.colorbar(image)colorbar.set_label('DTC (us/ft)', fontsize=14, fontweight='bold')plt.show()

Quando eseguiamo questo, otteniamo la seguente mappa che mostra la variazione della lentezza di compressione acustica per la Formazione Balder attraverso i nostri 65 pozzi.

Griglia di dati per la lentezza di compressione acustica (DTC) generata utilizzando pykrige. Immagine dell'autore.

Possiamo vedere che intorno ai 59-60 gradi N abbiamo rocce molto più veloci, e nelle regioni nord-est e sud-ovest abbiamo rocce molto più lente.

Per interpretarlo, dovremmo capire a quale profondità è la formazione in ciascuno di questi pozzi. Questo ci permetterà di identificare se la differenza è correlata alla sepoltura e alla compattazione o ad altri processi geologici.

Vedremo come fare questo in un futuro articolo.

Visualizzazione dell’incertezza del Kriging

Una delle cose chiave quando si guarda dati come questo è capire l’incertezza associata al kriging.

Possiamo farlo riutilizzando lo stesso codice di plotting, e invece di passare zstar, possiamo sostituirlo con la variabile ss che abbiamo creato in precedenza.

fig, ax = plt.subplots(figsize=(10,10))image = ax.imshow(ss, extent=(1.5, 4.5, 57.5, 62), origin='lower')ax.set_xlabel('Longitudine', fontsize=14, fontweight='bold')ax.set_ylabel('Latitudine', fontsize=14, fontweight='bold')scatter = ax.scatter(x=df['LON'], y=df['LAT'], color='black')colorbar = fig.colorbar(image)colorbar.set_label('DTC (us/ft)', fontsize=14, fontweight='bold')plt.show()

Con il seguente grafico, siamo in grado di vedere le aree in cui abbiamo un’alta o bassa incertezza.

Una griglia di dati di incertezza per la lentezza di compressione acustica (DTC) generata utilizzando pykrige. Immagine dell'autore.

In aree in cui abbiamo una copertura minore dai pozzi, avremo una maggiore incertezza, mentre in aree in cui abbiamo più pozzi, la nostra incertezza sarà molto inferiore.

Sommario

In questo tutorial, abbiamo visto come possiamo prendere i valori medi per una misura di registro di pozzo (DTC) e mapparli su un’intera regione. Ciò ci consente di comprendere le tendenze dei nostri dati su un’area geografica.

Tuttavia, quando guardiamo questi dati, dobbiamo tenere a mente che stiamo guardando una superficie 2D invece di una struttura 3D più complessa, che incontriamo all’interno del sottosuolo. Pertanto, le variazioni nella misurazione potrebbero essere attribuibili a variazioni nella profondità.

Dataset utilizzato

Il dataset utilizzato in questo articolo è un sottoinsieme di un dataset di formazione utilizzato come parte di una competizione di Machine Learning organizzata da Xeek e FORCE 2020 (Bormann et al., 2020). È rilasciato con licenza NOLD 2.0 dal governo norvegese, i dettagli della quale possono essere trovati qui: Licenza norvegese per dati governativi aperti (NLOD) 2.0. È possibile accedere al dataset completo qui.

Il riferimento completo per il dataset è:

Bormann, Peter, Aursand, Peder, Dilib, Fahad, Manral, Surrender, & Dischington, Peter. (2020). FORCE 2020 Well well log and lithofacies dataset for machine learning competition [Data set]. Zenodo. http://doi.org/10.5281/zenodo.4351156

Grazie per la lettura. Prima di andare, dovresti assolutamente iscriverti al mio contenuto e ricevere i miei articoli nella tua casella di posta. Puoi farlo qui!

In secondo luogo, puoi ottenere la piena esperienza Nisoo e supportare migliaia di altri scrittori e me iscrivendoti a una membership. Costa solo $ 5 al mese e hai pieno accesso a tutti i fantastici articoli di Nisoo, oltre alla possibilità di guadagnare denaro con la tua scrittura.

Se ti registri usando il mio link, mi supporterai direttamente con una parte della tua quota e non ti costerà di più. Se lo fai, grazie mille per il tuo supporto.