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

.. only:: html

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

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

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

.. _sphx_glr_auto_examples_example_dss_line.py:


Remove line noise with ZapLine
==============================

Find a spatial filter to get rid of line noise [1]_.

Uses meegkit.dss_line().

References
----------
.. [1] de Cheveigné, A. (2019). ZapLine: A simple and effective method to
    remove power line artifacts [Preprint]. https://doi.org/10.1101/782029

.. GENERATED FROM PYTHON SOURCE LINES 15-26

.. code-block:: Python

    # Authors: Maciej Szul <maciej.szul@isc.cnrs.fr>
    #          Nicolas Barascud <nicolas.barascud@gmail.com>
    import os

    import matplotlib.pyplot as plt
    import numpy as np
    from scipy import signal

    from meegkit import dss
    from meegkit.utils import create_line_data, unfold








.. GENERATED FROM PYTHON SOURCE LINES 27-29

Line noise removal
=============================================================================

.. GENERATED FROM PYTHON SOURCE LINES 31-34

Remove line noise with dss_line()
-----------------------------------------------------------------------------
We first generate some noisy data to work with

.. GENERATED FROM PYTHON SOURCE LINES 34-45

.. code-block:: Python

    sfreq = 250
    fline = 50
    nsamples = 10000
    nchans = 10
    data = create_line_data(n_samples=3 * nsamples, n_chans=nchans,
                            n_trials=1, fline=fline / sfreq, SNR=2)[0]
    data = data[..., 0]  # only take first trial

    # Apply dss_line (ZapLine)
    out, _ = dss.dss_line(data, fline, sfreq, nkeep=1)





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

 .. code-block:: none

    [PCA] Explained variance of selected components : 85.44%
    Power of components removed by DSS: 0.82




.. GENERATED FROM PYTHON SOURCE LINES 46-47

Plot before/after

.. GENERATED FROM PYTHON SOURCE LINES 47-60

.. code-block:: Python

    f, ax = plt.subplots(1, 2, sharey=True)
    f, Pxx = signal.welch(data, sfreq, nperseg=500, axis=0, return_onesided=True)
    ax[0].semilogy(f, Pxx)
    f, Pxx = signal.welch(out, sfreq, nperseg=500, axis=0, return_onesided=True)
    ax[1].semilogy(f, Pxx)
    ax[0].set_xlabel("frequency [Hz]")
    ax[1].set_xlabel("frequency [Hz]")
    ax[0].set_ylabel("PSD [V**2/Hz]")
    ax[0].set_title("before")
    ax[1].set_title("after")
    plt.show()





.. image-sg:: /auto_examples/images/sphx_glr_example_dss_line_001.png
   :alt: before, after
   :srcset: /auto_examples/images/sphx_glr_example_dss_line_001.png
   :class: sphx-glr-single-img





.. GENERATED FROM PYTHON SOURCE LINES 61-64

Remove line noise with dss_line_iter()
-----------------------------------------------------------------------------
We first load some noisy data to work with

.. GENERATED FROM PYTHON SOURCE LINES 64-72

.. code-block:: Python

    data = np.load(os.path.join("..", "tests", "data", "dss_line_data.npy"))
    fline = 50
    sfreq = 200
    print(data.shape)  # n_samples, n_chans, n_trials

    # Apply dss_line(), removing only one component
    out1, _ = dss.dss_line(data, fline, sfreq, nremove=1, nfft=400)





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

 .. code-block:: none

    (400, 64, 267)
    Power of components removed by DSS: 0.00




.. GENERATED FROM PYTHON SOURCE LINES 73-75

Now try dss_line_iter(). This applies dss_line() repeatedly until the
artifact is gone

.. GENERATED FROM PYTHON SOURCE LINES 75-78

.. code-block:: Python

    out2, iterations = dss.dss_line_iter(data, fline, sfreq, nfft=400)
    print(f"Removed {iterations} components")





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

 .. code-block:: none

    Power of components removed by DSS: 0.00
    Iteration 0 score: 8.822795627883339e-30
    Power of components removed by DSS: 0.00
    Iteration 1 score: 5.459919539000629e-30
    Power of components removed by DSS: 0.00
    Iteration 2 score: 1.8597287845371137e-31
    Power of components removed by DSS: 0.00
    Iteration 3 score: -2.8615484850989943e-31
    Removed 3 components




.. GENERATED FROM PYTHON SOURCE LINES 79-80

Plot results with dss_line() vs. dss_line_iter()

.. GENERATED FROM PYTHON SOURCE LINES 80-94

.. code-block:: Python

    f, ax = plt.subplots(1, 2, sharey=True)
    f, Pxx = signal.welch(unfold(out1), sfreq, nperseg=200, axis=0,
                          return_onesided=True)
    ax[0].semilogy(f, Pxx, lw=.5)
    f, Pxx = signal.welch(unfold(out2), sfreq, nperseg=200, axis=0,
                          return_onesided=True)
    ax[1].semilogy(f, Pxx, lw=.5)
    ax[0].set_xlabel("frequency [Hz]")
    ax[1].set_xlabel("frequency [Hz]")
    ax[0].set_ylabel("PSD [V**2/Hz]")
    ax[0].set_title("dss_line")
    ax[1].set_title("dss_line_iter")
    plt.tight_layout()
    plt.show()



.. image-sg:: /auto_examples/images/sphx_glr_example_dss_line_002.png
   :alt: dss_line, dss_line_iter
   :srcset: /auto_examples/images/sphx_glr_example_dss_line_002.png
   :class: sphx-glr-single-img






.. rst-class:: sphx-glr-timing

   **Total running time of the script:** (0 minutes 10.437 seconds)


.. _sphx_glr_download_auto_examples_example_dss_line.py:

.. only:: html

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

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

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

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

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


.. only:: html

 .. rst-class:: sphx-glr-signature

    `Gallery generated by Sphinx-Gallery <https://sphinx-gallery.github.io>`_