Merge pull request #44 from MiBiPreT/43-extend-examples

43 extend examples

examples/ex01_Griftpark/example_01_grift.py → examples/ex01_Griftpark/example_01_grift_NAscreening.py

@@ -70,4 +70,4 @@
 ### Create activity plot linking contaminant concentration to metabolite occurence
 ### and NA screening
 
-fig, ax = activity(data)
+fig, ax = activity(data,save_fig='grift_NA_activity.png',dpi = 300)

examples/ex01_Griftpark/example_01_grift_isotope.py

@@ -0,0 +1,88 @@
+#!/usr/bin/env python3
+# -*- coding: utf-8 -*-
+"""Lambda plot of isotopes.
+
+Script reproducing Figure 3.10 in the PhD thesis of Suzanne Faber, 2023
+'Field investigations and reactive transport modelling of biodegrading
+ coal tar compounds at a complex former manufactured gas plant''
+
+data provided on personal basis
+
+@author: Alraune Zech
+"""
+
+import matplotlib.pyplot as plt
+import numpy as np
+from mibipret.analysis.reduction.stable_isotope_regression import Lambda_regression
+from mibipret.analysis.reduction.stable_isotope_regression import extract_isotope_data
+from mibipret.data.data import load_csv
+from mibipret.data.data import standardize
+
+plt.close('all')
+
+# import sys
+# path = '/home/alraune/GitHub/MiBiPreT/mibipret/mibipret/'
+# sys.path.append(path)
+# from data.data import load_csv
+# from data.data import standardize
+# from analysis.reduction.stable_isotope_regression import Lambda_regression
+# from analysis.reduction.stable_isotope_regression import extract_isotope_data
+
+###------------------------------------------------------------------------###
+### Script settings
+verbose = True
+# molecules = ['benzene','toluene','ethylbenzene','pm_xylene','naphthalene','indane'] #make sure to use standard names
+molecules = ['Benzene','Toluene','Ethylbenzene','m,p-Xylene','Indane','Naphthalene']
+
+###------------------------------------------------------------------------###
+### File path settings
+file_csv = "./grift_BTEXIIN_isotopes.csv"
+save_fig = "./example_01_grift_isotopes.pdf"
+
+###------------------------------------------------------------------------###
+### Load and standardize data of isotopes
+
+isotopes_raw,units = load_csv(file_csv,
+                                  verbose = verbose)
+
+isotopes,units = standardize(isotopes_raw,
+                                reduce = True,
+                                verbose=verbose)
+
+###------------------------------------------------------------------------###
+### Figure settings
+
+wells = ['C','B2','B']
+mto = ['s','o','^','<','>','d','*']
+
+fig, axes = plt.subplots(figsize=[7.5,9],ncols = 2, nrows = 3)
+ax = axes.flat
+
+for j,molecule in enumerate(molecules):
+
+    for i,well in enumerate(wells):
+        data = isotopes[isotopes["obs_well"] == well]
+
+        x,y = extract_isotope_data(data,molecule)
+        results = Lambda_regression(x,y,validate_indices = True)
+
+        # Plot the scatter plot
+        ax[j].scatter(x, y, marker=mto[i],zorder = 3)
+        if len(x)>2:
+            # Create a trendline
+            polynomial = np.poly1d(results['coefficients'])
+            trendline_x = np.linspace(np.min(results['delta_C']), np.max(results['delta_C']), 50)
+            trendline_y = polynomial(trendline_x)
+            ax[j].plot(trendline_x, trendline_y,
+                       label=r'{}, $\Lambda = {:.0f}$'.format(well,results['coefficients'][0]))
+
+        ax[j].set_title(molecule)
+        ax[j].grid(True,zorder = 0)
+        ax[j].legend(title = r'Well & Lambda:')
+        if j%2 == 0:
+            ax[j].set_ylabel(r'$\delta^2$H')
+        if j >= len(ax)-2:
+            ax[j].set_xlabel(r'$\delta^{{13}}$C')
+
+fig.tight_layout()
+# plt.savefig(save_fig)

examples/ex01_Griftpark/grift_BTEXIIN_isotopes.csv

@@ -0,0 +1,17 @@
+sample_I,Well,Depth,delta_2H-Benzene,delta_2H-Toluene,delta_2H-Ethylbenzene,"delta_2H-m,p-Xylene",delta_2H-o-Xylene+Styrene,delta_2H-Indane,delta_2H-Naphthalene,delta_13C-Benzene,delta_13C-Toluene,delta_13C-Ethylbenzene,"delta_13C-m,p-Xylene",delta_13C-o-Xylene,delta_13C-Indane,delta_13C-Naphthalene,delta_13C-Indene
+,,m,mUr,mUr,mUr,mUr,mUr,mUr,mUr,mUr,mUr,mUr,mUr,mUr,mUr,mUr,mUr
+c-15,C,15.5,-99,170,-165,-131,-108,-153,-52,-26.7,-18.8,-25.2,-26.2,-25.6,-23.8,-27,-21.1
+c-19,C,19,-67,,-167,-122,-101,-155,-49,-25.8,-16.9,-24.2,-24.7,-24.6,-24.8,-26.6,
+c-26,C,26,-93,139,-155,-134,-96,-147,-55,-26,-21.6,-24.3,-25.7,-24.6,-25.2,-26.7,-24
+c-38,C,38.75,-72,-46,-155,-127,-100,-109,-47,-26,-25.3,-24.4,-25.5,-24.7,-25.2,-26.8,-25.9
+c-45,C,45.75,-91,13,-163,-138,-113,-159,-56,-25.5,-24.9,-24.3,-25.4,-25.2,-25.1,-26.8,-26.3
+b-12,B,12,-106,,80,62,70,-142,-43,-26.1,-22.8,-21.5,-17.5,-14,-24.7,-26.3,-23.2
+b-15,B,15.5,-110,,-154,-99,-43,-145,-54,-25.8,-17.8,-23.2,-24.5,-21.1,-24.7,-26.8,
+b-17,B,17,-118,78,-169,-115,-64,-159,-39,-24.1,-22.1,-23.4,-25.1,-22.1,-24.9,-26.7,
+b-19,B,19,-117,,-170,-49,9,-140,-42,-24.1,-20.2,-24.2,-22,-18,-24,-26.7,
+b-20,B,20.5,-90,,-173,-65,-22,-136,-34,-25.3,-18.7,-23.7,-23.1,-19.4,-24.6,-26.7,-23
+b-28,B,28,,,-162,-58,-97,,-28,-20.9,-23.8,-24.2,-21,-23.8,-26.9,-26.3,-23.7
+b2-7,B2,7,-77,,-179,-165,-138,-132,-55,-26.5,-21,-24.1,-26.2,-25.7,-27.6,-27,-25.7
+B2-10,B2,10.5,-75,,-161,-136,,-150,-60,-26.2,-20.6,-22.8,-23.7,-18.1,-26.9,-26.8,-25.4
+B2-14,B2,14,-51,38,-164,,-126,-106,-27,-25.2,-23.8,-24.1,-26,-25.1,-25.6,-27.2,-25.8
+B2-15,B2,15.5,-61,,2,-121,-11,-81,-60,-25.7,-18.1,-22.5,-25,-21.1,-23,-26.9,-24.4

examples/ex02_Amersfoort/example_02_amersfoort_isotopes_lambda.py

@@ -0,0 +1,109 @@
+#!/usr/bin/env python3
+# -*- coding: utf-8 -*-
+"""Lambda plot of isotopes.
+
+Script reproducing figure 8 in Paper van Leeuwen et al., 2022
+'Anaerobic degradation of benzene and other aromatic hydrocarbons in a
+tar-derived plume: Nitrate versus iron reducing conditions', J. of Cont. Hydrol
+
+data provided on personal basis
+Note: data does not match to data provided along manuscript, but is reconstructed
+from different data files
+
+@author: Alraune Zech
+"""
+
+import matplotlib.pyplot as plt
+import numpy as np
+from mibipret.analysis.reduction.stable_isotope_regression import Lambda_regression
+from mibipret.analysis.reduction.stable_isotope_regression import extract_isotope_data
+from mibipret.data.data import load_excel
+from mibipret.data.data import standardize
+from mibipret.visualize.stable_isotope_plots import Lambda_plot
+
+# import sys
+# path = '/home/alraune/GitHub/MiBiPreT/mibipret/mibipret/'
+# sys.path.append(path)
+# from analysis.reduction.stable_isotope_regression import Lambda_regression
+# from analysis.reduction.stable_isotope_regression import extract_isotope_data
+# from visualize.stable_isotope_plots import Lambda_plot
+
+plt.close('all')
+
+###------------------------------------------------------------------------###
+### Script settings
+verbose = True
+molecules = ['benzene','toluene','ethylbenzene','pm_xylene','naphthalene','indene']
+molecules_analysis = []
+###------------------------------------------------------------------------###
+
+### File path settings
+file_path = './amersfoort.xlsx'
+# file_path = "./amersfoort_isotopes_match.csv"
+
+###------------------------------------------------------------------------###
+### Load and standardize data of isotopes
+# isotopes_raw,units = md.load_csv(file_path,
+#                                  verbose = verbose)
+
+isotopes_raw,units = load_excel(file_path,
+                                   sheet_name = 'isotopes',
+                                   verbose = verbose)
+
+isotopes,units = standardize(isotopes_raw,
+                                reduce = True,
+                                verbose=verbose)
+
+###------------------------------------------------------------------------###
+### Lambda regression and Lambda regression plot for all molecules (separately)
+
+for j,molecule in enumerate(molecules):
+
+    C_data,H_data = extract_isotope_data(isotopes,molecule)
+
+    results = Lambda_regression(C_data,
+                                H_data,
+                                validate_indices = True,
+                                verbose = verbose,
+                                )
+
+    molecules_analysis.append(results)
+
+    Lambda_plot(**results)#,save_fig = 'Amersfoort_isotope_Lambda_{}.pdf'.format(molecule))
+    # Lambda_plot(**results,save_fig = 'Amersfoort_isotope_Lambda_{}.png'.format(molecule),dpi = 300)
+
+###------------------------------------------------------------------------###
+### Lambda regression and plot of results of all molecules in one plot
+
+plt.figure(j)
+fig, axes = plt.subplots(figsize=[7.5,9],
+                          ncols = 2,
+                          nrows = 3)
+ax = axes.flat
+for j,molecule in enumerate(molecules):
+
+    results =  molecules_analysis[j]
+
+    x,y = results['delta_C'], results['delta_H']
+
+    # Plot the scatter plot
+    ax[j].scatter(x,y, zorder = 3,label = 'data')
+    ax[j].set_title(molecule)
+    ax[j].grid(True,zorder = 0)
+
+    # plot trendlines
+    polynomial = np.poly1d(results['coefficients'])
+    trendline_x = np.linspace(np.min(x), np.max(x), 50)
+    trendline_y = polynomial(trendline_x)
+    ax[j].plot(trendline_x, trendline_y, label=r'$\Lambda = {:.0f}$'.format(results['coefficients'][0]))
+
+    if j%2 == 0:
+        ax[j].set_ylabel(r'$\delta^2$H')
+
+    if j >= len(ax)-2:
+        ax[j].set_xlabel(r'$\delta^{{13}}$C')
+
+    ax[j].legend()
+
+fig.tight_layout()
+# plt.savefig("Amersfoort_isotope_Lambda.pdf")

mibipret/analysis/reduction/stable_isotope_regression.py

@@ -5,6 +5,8 @@
 @author: Alraune Zech
 """
 import numpy as np
+from mibipret.data.names import names_contaminants
+from mibipret.data.names import names_isotopes
 
 
 def Lambda_regression(delta_C,
@@ -343,3 +345,53 @@ def valid_indices(data1,
         valid_indices *= (data1 != 0) & (data2 != 0)
 
     return data1[valid_indices],data2[valid_indices]
+
+def extract_isotope_data(df,
+                         molecule,
+                         name_13C = 'delta_13C',
+                         name_2H = 'delta_2H',
+                         ):
+    """Extracts isotope data from standardised input-dataframe.
+
+    Parameters
+    ----------
+    df : pd.dataframe
+        numeric (observational) data
+    molecule : str
+        name of contaminant molecule to extract isotope data for
+    name_13C : str, default 'delta_13C' (standard name)
+        name of C13 isotope to extract data for
+    name_2H : str, default 'delta_2H' (standard name)
+        name of deuterium isotope to extract data for
+
+    Returns
+    -------
+    C_data : np.array
+        numeric isotope data
+    H_data : np.array
+        numeric isotope data
+
+    """
+    molecule_standard = names_contaminants.get(molecule.lower(), False)
+    isotope_13C = names_isotopes.get(name_13C.lower(), False)
+    isotope_2H = names_isotopes.get(name_2H.lower(), False)
+
+    if molecule_standard is False:
+        raise ValueError("Contaminant (name) unknown: {}".format(molecule))
+    if isotope_13C is False:
+        raise ValueError("Isotope (name) unknown: {}".format(name_13C))
+    if isotope_2H is False:
+        raise ValueError("Isotope (name) unknown: {}".format(name_2H))
+
+    name_C = '{}-{}'.format(isotope_13C,molecule_standard)
+    name_H = '{}-{}'.format(isotope_2H,molecule_standard)
+
+    if name_C not in df.columns.to_list():
+        raise ValueError("No isotope data available for : {}".format(name_C))
+    if name_H not in df.columns.to_list():
+        raise ValueError("No isotope data available for : {}".format(name_H))
+
+    C_data = df[name_C].values
+    H_data = df[name_H].values
+
+    return C_data, H_data

mibipret/visualize/activity.py

@@ -26,6 +26,7 @@
     ec = 'k',
     lw = 0.5,
     loc = 'lower right',
+    dpi = 300,
     )
 
 def activity(
@@ -149,7 +150,7 @@ def activity(
     ### Save figure to file if file path provided
     if save_fig is not False:
         try:
-            plt.savefig(save_fig)
+            plt.savefig(save_fig,dpi = settings['dpi'])
             print("Save Figure to file:\n", save_fig)
         except OSError:
             print("WARNING: Figure could not be saved. Check provided file path and name: {}".format(save_fig))

mibipret/visualize/ordination_plot.py

@@ -11,7 +11,6 @@
 
 DEF_settings = dict(
     figsize = [3.75,3.75],
-    dpi = 300,
     label_fontsize = 8,
     loading_fontsize = 8,
     score_fontsize = 6,
@@ -29,6 +28,7 @@
     score_facecolor = 'none',
     score_marker = 'o',
     score_marker_size = 45,
+    dpi = 300,
     )
 
 def ordination_plot(ordination_output,
@@ -151,9 +151,7 @@ def ordination_plot(ordination_output,
 
     ### ---------------------------------------------------------------------------
     ### Create Figure, finally!
-    fig, ax = plt.subplots(figsize=settings['figsize'],
-                           dpi=settings['dpi'],
-                           )
+    fig, ax = plt.subplots(figsize=settings['figsize'])
     texts = []
 
     # Plotting the ordination scores by iterating over every coordinate
@@ -245,7 +243,7 @@ def ordination_plot(ordination_output,
     plt.tight_layout()
     if save_fig is not False:
         try:
-            plt.savefig(save_fig)
+            plt.savefig(save_fig,dpi = settings['dpi'])
             print("Figure saved to file:\n", save_fig)
         except OSError:
             print("WARNING: Figure could not be saved. Check provided file path and name: {}".format(save_fig))

mibipret/visualize/stable_isotope_plots.py

@@ -17,6 +17,7 @@
     fit_color = 'red',
     intercept_color = 'forestgreen',
     loc = 'best',
+    dpi = 300,
     )
 
 def Lambda_plot(delta_C,
@@ -115,7 +116,7 @@ def Lambda_plot(delta_C,
     ### Save figure to file if file path provided
     if save_fig is not False:
         try:
-            plt.savefig(save_fig)
+            plt.savefig(save_fig,dpi = settings['dpi'])
             print("Save Figure to file:\n", save_fig)
         except OSError:
             print("WARNING: Figure could not be saved. Check provided file path and name: {}".format(save_fig))
@@ -207,7 +208,7 @@ def Rayleigh_fractionation_plot(concentration,
     ### Save figure to file if file path provided
     if save_fig is not False:
         try:
-            plt.savefig(save_fig)
+            plt.savefig(save_fig,dpi = settings['dpi'])
             print("Save Figure to file:\n", save_fig)
         except OSError:
             print("WARNING: Figure could not be saved. Check provided file path and name: {}".format(save_fig))
@@ -314,7 +315,7 @@ def Keeling_plot(concentration,
     ### Save figure to file if file path provided
     if save_fig is not False:
         try:
-            plt.savefig(save_fig)
+            plt.savefig(save_fig,dpi = settings['dpi'])
             print("Save Figure to file:\n", save_fig)
         except OSError:
             print("WARNING: Figure could not be saved. Check provided file path and name: {}".format(save_fig))

tests/test_data.py

@@ -15,6 +15,7 @@
 from mibipret.data.data import standardize
 
 path_data = "./mibipret/data"
+#path_data = "../mibipret/data"
 
 class TestData:
     """Class for testing data module of mibipret."""