.. DO NOT EDIT.
.. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY.
.. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE:
.. "auto_data_analysis/estimate_stochastic_processes/plot_estimate_non_stationary_covariance_model.py"
.. LINE NUMBERS ARE GIVEN BELOW.

.. only:: html

    .. note::
        :class: sphx-glr-download-link-note

        :ref:`Go to the end <sphx_glr_download_auto_data_analysis_estimate_stochastic_processes_plot_estimate_non_stationary_covariance_model.py>`
        to download the full example code.

.. rst-class:: sphx-glr-example-title

.. _sphx_glr_auto_data_analysis_estimate_stochastic_processes_plot_estimate_non_stationary_covariance_model.py:


Estimate a non stationary covariance function
=============================================

.. GENERATED FROM PYTHON SOURCE LINES 7-29

The objective of this use case is to estimate :math:`C` from several
fields generated by the process :math:`X`. We suppose that the process
is not stationary.

In the following example, we illustrate the estimation of the non stationary covariance model

:math:`C : \mathcal{D} \times [-4, 4] \rightarrow \mathbb{R}` defined by:

.. math::
   \begin{aligned}
     \displaystyle C(s,t) = \exp\left(-\dfrac{4|s-t|}{1+s^2+t^2}\right)
   \end{aligned}

The domain :math:`\mathcal{D}` is discretized on a mesh :math:`\mathcal{M}` which is a time grid with 64 points.
We build a Normal process :math:`X: \Omega \times [-4, 4]  \rightarrow \mathbb{R}` with zero mean and
:math:`C` as covariance function.
We discretize the covariance model :math:`C` using :math:`C(t_k, t_\ell)`
for each :math:`(t_k, t_\ell)\in \mathcal{M} \times \mathcal{M}`.
We get a :math:`N=10^3` fields from the process :math:`X` from which we
estimate the covariance model :math:`C`.

We use the object :class:`~openturns.NonStationaryCovarianceModelFactory` which creates a :class:`~openturns.UserDefinedCovarianceModel`.

.. GENERATED FROM PYTHON SOURCE LINES 31-37

.. code-block:: Python

    import math as m
    import openturns as ot
    import openturns.viewer as viewer
    from matplotlib import pyplot as plt









.. GENERATED FROM PYTHON SOURCE LINES 38-39

Create the time grid

.. GENERATED FROM PYTHON SOURCE LINES 39-46

.. code-block:: Python

    t0 = -4.0
    tmax = 4.0
    N = 64
    dt = (tmax - t0) / N
    tgrid = ot.RegularGrid(t0, dt, N)









.. GENERATED FROM PYTHON SOURCE LINES 47-48

Create the covariance function at (s,t)

.. GENERATED FROM PYTHON SOURCE LINES 48-52

.. code-block:: Python

    def C(s, t):
        return m.exp(-4.0 * abs(s - t) / (1 + (s * s + t * t)))









.. GENERATED FROM PYTHON SOURCE LINES 53-54

Draw...

.. GENERATED FROM PYTHON SOURCE LINES 54-66

.. code-block:: Python

    def f(X):
        s, t = X
        return [C(s, t)]


    func = ot.PythonFunction(2, 1, f)
    func.setDescription([":math:`s`", ":math:`t`", ":math:`cov`"])
    graph = func.draw([t0] * 2, [tmax] * 2)
    graph.setTitle("Original covariance model")
    graph.setLegendPosition("")
    view = viewer.View(graph)




.. image-sg:: /auto_data_analysis/estimate_stochastic_processes/images/sphx_glr_plot_estimate_non_stationary_covariance_model_001.svg
   :alt: Original covariance model
   :srcset: /auto_data_analysis/estimate_stochastic_processes/images/sphx_glr_plot_estimate_non_stationary_covariance_model_001.svg
   :class: sphx-glr-single-img





.. GENERATED FROM PYTHON SOURCE LINES 67-68

Create data from a non stationary Normal process Omega * [0,T]--> R

.. GENERATED FROM PYTHON SOURCE LINES 68-79

.. code-block:: Python


    # Create the collection of HermitianMatrix
    covariance = ot.CovarianceMatrix(N)
    for k in range(N):
        s = tgrid.getValue(k)
        for ll in range(k + 1):
            t = tgrid.getValue(ll)
            covariance[k, ll] = C(s, t)

    covmodel = ot.UserDefinedCovarianceModel(tgrid, covariance)








.. GENERATED FROM PYTHON SOURCE LINES 80-82

Create the Normal process with that covariance model
based on the mesh tgrid

.. GENERATED FROM PYTHON SOURCE LINES 82-88

.. code-block:: Python

    process = ot.GaussianProcess(covmodel, tgrid)

    # Get a sample of fields from the process
    N = 1000
    sample = process.getSample(N)








.. GENERATED FROM PYTHON SOURCE LINES 89-90

The covariance model factory

.. GENERATED FROM PYTHON SOURCE LINES 90-95

.. code-block:: Python

    factory = ot.NonStationaryCovarianceModelFactory()

    # Estimation on a sample
    estimatedModel = factory.build(sample)








.. GENERATED FROM PYTHON SOURCE LINES 96-101

.. code-block:: Python

    graph = estimatedModel.draw(0, 0, t0, tmax, 256, False)
    graph.setTitle("Estimated covariance model")
    graph.setLegendPosition("")
    view = viewer.View(graph)
    plt.show()



.. image-sg:: /auto_data_analysis/estimate_stochastic_processes/images/sphx_glr_plot_estimate_non_stationary_covariance_model_002.svg
   :alt: Estimated covariance model
   :srcset: /auto_data_analysis/estimate_stochastic_processes/images/sphx_glr_plot_estimate_non_stationary_covariance_model_002.svg
   :class: sphx-glr-single-img






.. _sphx_glr_download_auto_data_analysis_estimate_stochastic_processes_plot_estimate_non_stationary_covariance_model.py:

.. only:: html

  .. container:: sphx-glr-footer sphx-glr-footer-example

    .. container:: sphx-glr-download sphx-glr-download-jupyter

      :download:`Download Jupyter notebook: plot_estimate_non_stationary_covariance_model.ipynb <plot_estimate_non_stationary_covariance_model.ipynb>`

    .. container:: sphx-glr-download sphx-glr-download-python

      :download:`Download Python source code: plot_estimate_non_stationary_covariance_model.py <plot_estimate_non_stationary_covariance_model.py>`

    .. container:: sphx-glr-download sphx-glr-download-zip

      :download:`Download zipped: plot_estimate_non_stationary_covariance_model.zip <plot_estimate_non_stationary_covariance_model.zip>`