.. 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_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_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_stationary_covariance_model.py:


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

.. GENERATED FROM PYTHON SOURCE LINES 7-21

The objective here is to estimate a stationary covariance model from data.

The library builds an estimation of the stationary covariance function on a
*ProcessSample* or *TimeSeries* using the previous algorithm implemented in
the *StationaryCovarianceModelFactory* class.
The result consists in a *UserDefinedStationaryCovarianceModel* which is easy to manipulate.

Such an object is composed of a time grid and a collection of :math:`K`
square matrices of dimension d.
:math:`K` corresponds to the number of time steps of the final time grid on
which the covariance is estimated.
When estimated from a time series , the *UserDefinedStationaryCovarianceModel*
may have a time grid different from the initial time grid of the time series.


.. GENERATED FROM PYTHON SOURCE LINES 23-29

.. code-block:: Python

    import openturns as ot
    import openturns.viewer as viewer
    from matplotlib import pylab as plt

    ot.Log.Show(ot.Log.NONE)








.. GENERATED FROM PYTHON SOURCE LINES 30-31

Create some 1-d normal process data with an Exponential covariance model

.. GENERATED FROM PYTHON SOURCE LINES 31-54

.. code-block:: Python


    # Dimension parameter
    dim = 1

    # Create the time grid
    t0 = 0.0
    N = 300
    t1 = 20.0
    dt = (t1 - t0) / N
    tgrid = ot.RegularGrid(t0, dt, N)

    # Create the covariance model
    amplitude = [1.0] * dim
    scale = [1.0] * dim
    covmodel = ot.ExponentialModel(scale, amplitude)

    # Create a stationary Normal process with that covariance model
    process = ot.GaussianProcess(covmodel, tgrid)

    # Create a time series and a sample of time series
    tseries = process.getRealization()
    sample = process.getSample(1000)








.. GENERATED FROM PYTHON SOURCE LINES 55-56

Build a factory of stationary covariance function

.. GENERATED FROM PYTHON SOURCE LINES 56-66

.. code-block:: Python

    covarianceFactory = ot.StationaryCovarianceModelFactory()

    # Set the spectral factory algorithm
    segmentNumber = 5
    spectralFactory = ot.WelchFactory(ot.Hann(), segmentNumber)
    covarianceFactory.setSpectralModelFactory(spectralFactory)

    # Check the current spectral factory
    print(covarianceFactory.getSpectralModelFactory())





.. rst-class:: sphx-glr-script-out

 .. code-block:: none

    class=WelchFactory window = class=FilteringWindows implementation=class=Hann blockNumber = 5 overlap = 0.5




.. GENERATED FROM PYTHON SOURCE LINES 67-68

Case 1 :  Estimation on a ProcessSample

.. GENERATED FROM PYTHON SOURCE LINES 68-85

.. code-block:: Python


    # The spectral model factory computes the spectral density function
    # without using the block and overlap arguments of the Welch factories
    estimatedModel_PS = covarianceFactory.build(sample)

    # Case 2 :  Estimation on a TimeSeries

    # The spectral model factory compute the spectral density function using
    # the block and overlap arguments of spectral model factories
    estimatedModel_TS = covarianceFactory.build(tseries)

    # Evaluate the covariance function at each time step
    # Care : if estimated from a time series, the time grid has changed
    for i in range(N):
        tau = tgrid.getValue(i)
        cov = estimatedModel_PS(tau)








.. GENERATED FROM PYTHON SOURCE LINES 86-87

Drawing...

.. GENERATED FROM PYTHON SOURCE LINES 87-108

.. code-block:: Python

    sampleValueEstimated = ot.Sample(N, 1)
    sampleValueModel = ot.Sample(N, 1)
    for i in range(N):
        t = tgrid.getValue(i)
        for j in range(i - 1):
            s = tgrid.getValue(j)
            estimatedValue = estimatedModel_PS(t, s)
            modelValue = covmodel(t, s)
            if j == 0:
                sampleValueEstimated[i, 0] = estimatedValue[0, 0]
                sampleValueModel[i, 0] = modelValue[0, 0]
    sampleT = tgrid.getVertices()
    graph = ot.Graph("Covariance estimation", "time", "Covariance value C(0,t)", True)
    curveEstimated = ot.Curve(sampleT, sampleValueEstimated, "Estimated model")
    graph.add(curveEstimated)
    curveModel = ot.Curve(sampleT, sampleValueModel, "Exact model")
    curveModel.setColor("red")
    graph.add(curveModel)
    graph.setLegendPosition("upper right")
    view = viewer.View(graph)
    plt.show()



.. image-sg:: /auto_data_analysis/estimate_stochastic_processes/images/sphx_glr_plot_estimate_stationary_covariance_model_001.png
   :alt: Covariance estimation
   :srcset: /auto_data_analysis/estimate_stochastic_processes/images/sphx_glr_plot_estimate_stationary_covariance_model_001.png
   :class: sphx-glr-single-img






.. _sphx_glr_download_auto_data_analysis_estimate_stochastic_processes_plot_estimate_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_stationary_covariance_model.ipynb <plot_estimate_stationary_covariance_model.ipynb>`

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

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