.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "examples/Models/plot_model_1d.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note :ref:`Go to the end ` to download the full example code. .. rst-class:: sphx-glr-example-title .. _sphx_glr_examples_Models_plot_model_1d.py: 1D Model with an infinite halfspace ----------------------------------- .. GENERATED FROM PYTHON SOURCE LINES 7-15 .. code-block:: Python from copy import deepcopy from geobipy import StatArray from geobipy import RectilinearMesh1D from geobipy import Model from geobipy import Distribution import matplotlib.pyplot as plt import numpy as np .. GENERATED FROM PYTHON SOURCE LINES 16-18 Instantiate the 1D Model with a Half Space ++++++++++++++++++++++++++++++++++++++++++ .. GENERATED FROM PYTHON SOURCE LINES 18-31 .. code-block:: Python # Make a test model with 10 layers, and increasing parameter values nLayers = 2 par = StatArray(np.linspace(0.001, 0.02, nLayers), "Conductivity", "$\\frac{S}{m}$") thk = StatArray(np.full(nLayers, fill_value=10.0)) thk[-1] = np.inf mesh = RectilinearMesh1D(widths = thk) mod = Model(mesh = mesh, values=par) plt.figure() mod.plot_grid(transpose=True, flip=True) .. image-sg:: /examples/Models/images/sphx_glr_plot_model_1d_001.png :alt: plot model 1d :srcset: /examples/Models/images/sphx_glr_plot_model_1d_001.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 32-39 Randomness and Model Perturbations ++++++++++++++++++++++++++++++++++ We can set the priors on the 1D model by assigning minimum and maximum layer depths and a maximum number of layers. These are used to create priors on the number of cells in the model, a new depth interface, new parameter values and the vertical gradient of those parameters. The halfSpaceValue is used as a reference value for the parameter prior. .. GENERATED FROM PYTHON SOURCE LINES 39-53 .. code-block:: Python from numpy.random import Generator from numpy.random import PCG64DXSM generator = PCG64DXSM(seed=0) prng = Generator(generator) # Set the priors mod.set_priors(value_mean=0.01, min_edge=1.0, max_edge=150.0, max_cells=30, solve_value=True, solve_gradient=True, prng=prng) .. GENERATED FROM PYTHON SOURCE LINES 54-55 We can evaluate the prior of the model using depths only .. GENERATED FROM PYTHON SOURCE LINES 55-59 .. code-block:: Python print('Log probability of the Model given its priors: ', mod.probability(False, False)) # Or with priors on its parameters, and parameter gradient with depth. print('Log probability of the Model given its priors: ', mod.probability(True, True)) .. rst-class:: sphx-glr-script-out .. code-block:: none Log probability of the Model given its priors: -3.367295829986474 Log probability of the Model given its priors: -9.157731937350919 .. GENERATED FROM PYTHON SOURCE LINES 60-62 To propose new models, we specify the probabilities of creating, removing, perturbing, and not changing a layer interface .. GENERATED FROM PYTHON SOURCE LINES 62-65 .. code-block:: Python pProposal = Distribution('LogNormal', 0.01, np.log(2.0)**2.0, linearSpace=True, prng=prng) mod.set_proposals(probabilities=[0.25, 0.25, 0.5, 0.25], proposal=pProposal, prng=prng) .. GENERATED FROM PYTHON SOURCE LINES 66-67 We can then perturb the layers of the model .. GENERATED FROM PYTHON SOURCE LINES 67-69 .. code-block:: Python remapped, perturbed = mod.perturb() .. GENERATED FROM PYTHON SOURCE LINES 70-76 .. code-block:: Python fig = plt.figure(figsize=(8, 6)) ax = plt.subplot(121) mod.pcolor(transpose=True, flip=True, log=10) # , grid=True) ax = plt.subplot(122) perturbed.pcolor(transpose=True, flip=True, log=10) # , grid=True) .. image-sg:: /examples/Models/images/sphx_glr_plot_model_1d_002.png :alt: plot model 1d :srcset: /examples/Models/images/sphx_glr_plot_model_1d_002.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none .. GENERATED FROM PYTHON SOURCE LINES 77-78 We can evaluate the prior of the model using depths only .. GENERATED FROM PYTHON SOURCE LINES 78-83 .. code-block:: Python print('Log probability of the Model given its priors: ',perturbed.probability(False, False)) # Or with priors on its parameters, and parameter gradient with depth. print('Log probability of the Model given its priors: ',perturbed.probability(True, True)) .. rst-class:: sphx-glr-script-out .. code-block:: none Log probability of the Model given its priors: -3.367295829986474 Log probability of the Model given its priors: -8.559817917094882 .. GENERATED FROM PYTHON SOURCE LINES 84-98 Perturbing a model multiple times +++++++++++++++++++++++++++++++++ In the stochasitic inference process, we perturb the model structure, and parameter values, multiple times. Each time the model is perturbed, we can record its state in a posterior distribution. For a 1D model, the parameter posterior is a 2D hitmap with depth in one dimension and the parameter value in the other. We also attach a 1D histogram for the number of layers, and a 1D histogram for the locations of interfaces. Since we have already set the priors on the Model, we can set the posteriors based on bins from from the priors. .. GENERATED FROM PYTHON SOURCE LINES 98-103 .. code-block:: Python mod.set_posteriors() mod0 = deepcopy(mod) .. GENERATED FROM PYTHON SOURCE LINES 104-105 Now we randomly perturb the model, and update its posteriors. .. GENERATED FROM PYTHON SOURCE LINES 105-114 .. code-block:: Python mod.update_posteriors() for i in range(1001): remapped, perturbed = mod.perturb() # And update the model posteriors perturbed.update_posteriors() mod = perturbed .. GENERATED FROM PYTHON SOURCE LINES 115-120 We can now plot the posteriors of the model. Remember in this case, we are simply perturbing the model structure and parameter values The proposal for the parameter values is fixed and centred around a single value. fig = plt.figure(figsize=(8, 6)) .. GENERATED FROM PYTHON SOURCE LINES 120-161 .. code-block:: Python # plt.subplot(131) # mod.nCells.posterior.plot() # ax = plt.subplot(132) # mod.values.posterior.pcolor(cmap='gray_r', colorbar=False, flipY=True, logX=10) # plt.subplot(133, sharey=ax) # mod.mesh.edges.posterior.plot(transpose=True, flipY=True) # plt.figure() # mod.plot_posteriors(**{"cmap": 'gray_r', # "xscale": 'log', # "noColorbar": True, # "flipY": True, # 'credible_interval_kwargs':{'axis': 1, # 'reciprocate': True, # 'xscale': 'log'}}) # mod.par.posterior.plotCredibleIntervals(xscale='log', axis=1) fig = plt.figure(figsize=(8, 6)) # gs = fig.add_gridspec(nrows=1, ncols=1) mod.plot_posteriors(axes=fig, edges_kwargs = { "transpose":True, "flipY":True }, parameter_kwargs = { "cmap": 'gray_r', "xscale": 'log', "colorbar": False, "flipY": True, 'credible_interval_kwargs':{ 'reciprocate':True, # 'axis': 1, 'xscale': 'log' } }, best = mod) plt.show() .. image-sg:: /examples/Models/images/sphx_glr_plot_model_1d_003.png :alt: plot model 1d :srcset: /examples/Models/images/sphx_glr_plot_model_1d_003.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 21.040 seconds) .. _sphx_glr_download_examples_Models_plot_model_1d.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: plot_model_1d.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: plot_model_1d.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: plot_model_1d.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_