{
  "nbformat_minor": 0, 
  "nbformat": 4, 
  "cells": [
    {
      "execution_count": null, 
      "cell_type": "code", 
      "source": [
        "%matplotlib inline"
      ], 
      "outputs": [], 
      "metadata": {
        "collapsed": false
      }
    }, 
    {
      "source": [
        "\n\n# Background information on filtering\n\n\nHere we give some background information on filtering in general,\nand how it is done in MNE-Python in particular.\nRecommended reading for practical applications of digital\nfilter design can be found in Parks & Burrus [1]_ and\nIfeachor and Jervis [2]_, and for filtering in an\nM/EEG context we recommend reading Widmann *et al.* 2015 [7]_.\nTo see how to use the default filters in MNE-Python on actual data, see\nthe `tut_artifacts_filter` tutorial.\n\nProblem statement\n=================\n\nThe practical issues with filtering electrophysiological data are covered\nwell by Widmann *et al.* in [7]_, in a follow-up to an article where they\nconclude with this statement:\n\n    Filtering can result in considerable distortions of the time course\n    (and amplitude) of a signal as demonstrated by VanRullen (2011) [[3]_].\n    Thus, filtering should not be used lightly. However, if effects of\n    filtering are cautiously considered and filter artifacts are minimized,\n    a valid interpretation of the temporal dynamics of filtered\n    electrophysiological data is possible and signals missed otherwise\n    can be detected with filtering.\n\nIn other words, filtering can increase SNR, but if it is not used carefully,\nit can distort data. Here we hope to cover some filtering basics so\nusers can better understand filtering tradeoffs, and why MNE-Python has\nchosen particular defaults.\n\n\nFiltering basics\n================\n\nLet's get some of the basic math down. In the frequency domain, digital\nfilters have a transfer function that is given by:\n\n\\begin{align}H(z) &= \\frac{b_0 + b_1 z^{-1} + b_2 z^{-2} + ... + b_M z^{-M}}\n                 {1 + a_1 z^{-1} + a_2 z^{-2} + ... + a_N z^{-M}} \\\\\n         &= \\frac{\\sum_0^Mb_kz^{-k}}{\\sum_1^Na_kz^{-k}}\\end{align}\n\nIn the time domain, the numerator coefficients $b_k$ and denominator\ncoefficients $a_k$ can be used to obtain our output data\n$y(n)$ in terms of our input data $x(n)$ as:\n\n\\begin{align}:label: summations\n\n    y(n) &= b_0 x(n) + b_1 x(n-1) + ... + b_M x(n-M)\n            - a_1 y(n-1) - a_2 y(n - 2) - ... - a_N y(n - N)\\\\\n         &= \\sum_0^M b_k x(n-k) - \\sum_1^N a_k y(n-k)\\end{align}\n\nIn other words, the output at time $n$ is determined by a sum over:\n\n    1. The numerator coefficients $b_k$, which get multiplied by\n       the previous input $x(n-k)$ values, and\n    2. The denominator coefficients $a_k$, which get multiplied by\n       the previous output $y(n-k)$ values.\n\nNote that these summations in :eq:`summations` correspond nicely to\n(1) a weighted `moving average`_ and (2) an autoregression_.\n\nFilters are broken into two classes: FIR_ (finite impulse response) and\nIIR_ (infinite impulse response) based on these coefficients.\nFIR filters use a finite number of numerator\ncoefficients $b_k$ ($\\forall k, a_k=0$), and thus each output\nvalue of $y(n)$ depends only on the $M$ previous input values.\nIIR filters depend on the previous input and output values, and thus can have\neffectively infinite impulse responses.\n\nAs outlined in [1]_, FIR and IIR have different tradeoffs:\n\n    * A causal FIR filter can be linear-phase -- i.e., the same time delay\n      across all frequencies -- whereas a causal IIR filter cannot. The phase\n      and group delay characteristics are also usually better for FIR filters.\n    * IIR filters can generally have a steeper cutoff than an FIR filter of\n      equivalent order.\n    * IIR filters are generally less numerically stable, in part due to\n      accumulating error (due to its recursive calculations).\n\nIn MNE-Python we default to using FIR filtering. As noted in Widmann *et al.*\n2015 [7]_:\n\n    Despite IIR filters often being considered as computationally more\n    efficient, they are recommended only when high throughput and sharp\n    cutoffs are required (Ifeachor and Jervis, 2002 [2]_, p. 321),\n    ...FIR filters are easier to control, are always stable, have a\n    well-defined passband, can be corrected to zero-phase without\n    additional computations, and can be converted to minimum-phase.\n    We therefore recommend FIR filters for most purposes in\n    electrophysiological data analysis.\n\nWhen designing a filter (FIR or IIR), there are always tradeoffs that\nneed to be considered, including but not limited to:\n\n    1. Ripple in the pass-band\n    2. Attenuation of the stop-band\n    3. Steepness of roll-off\n    4. Filter order (i.e., length for FIR filters)\n    5. Time-domain ringing\n\nIn general, the sharper something is in frequency, the broader it is in time,\nand vice-versa. This is a fundamental time-frequency tradeoff, and it will\nshow up below.\n\nFIR Filters\n===========\n\nFirst we will focus first on FIR filters, which are the default filters used by\nMNE-Python.\n\n"
      ], 
      "cell_type": "markdown", 
      "metadata": {}
    }, 
    {
      "source": [
        "Designing FIR filters\n---------------------\nHere we'll try designing a low-pass filter, and look at trade-offs in terms\nof time- and frequency-domain filter characteristics. Later, in\n`tut_effect_on_signals`, we'll look at how such filters can affect\nsignals when they are used.\n\nFirst let's import some useful tools for filtering, and set some default\nvalues for our data that are reasonable for M/EEG data.\n\n"
      ], 
      "cell_type": "markdown", 
      "metadata": {}
    }, 
    {
      "execution_count": null, 
      "cell_type": "code", 
      "source": [
        "import numpy as np\nfrom scipy import signal, fftpack\nimport matplotlib.pyplot as plt\n\nfrom mne.time_frequency.tfr import morlet\nfrom mne.viz import plot_filter, plot_ideal_filter\n\nimport mne\n\nsfreq = 1000.\nf_p = 40.\nflim = (1., sfreq / 2.)  # limits for plotting"
      ], 
      "outputs": [], 
      "metadata": {
        "collapsed": false
      }
    }, 
    {
      "source": [
        "Take for example an ideal low-pass filter, which would give a value of 1 in\nthe pass-band (up to frequency $f_p$) and a value of 0 in the stop-band\n(down to frequency $f_s$) such that $f_p=f_s=40$ Hz here\n(shown to a lower limit of -60 dB for simplicity):\n\n"
      ], 
      "cell_type": "markdown", 
      "metadata": {}
    }, 
    {
      "execution_count": null, 
      "cell_type": "code", 
      "source": [
        "nyq = sfreq / 2.  # the Nyquist frequency is half our sample rate\nfreq = [0, f_p, f_p, nyq]\ngain = [1, 1, 0, 0]\n\nthird_height = np.array(plt.rcParams['figure.figsize']) * [1, 1. / 3.]\nax = plt.subplots(1, figsize=third_height)[1]\nplot_ideal_filter(freq, gain, ax, title='Ideal %s Hz lowpass' % f_p, flim=flim)"
      ], 
      "outputs": [], 
      "metadata": {
        "collapsed": false
      }
    }, 
    {
      "source": [
        "This filter hypothetically achieves zero ripple in the frequency domain,\nperfect attenuation, and perfect steepness. However, due to the discontunity\nin the frequency response, the filter would require infinite ringing in the\ntime domain (i.e., infinite order) to be realized. Another way to think of\nthis is that a rectangular window in frequency is actually sinc_ function\nin time, which requires an infinite number of samples, and thus infinite\ntime, to represent. So although this filter has ideal frequency suppression,\nit has poor time-domain characteristics.\n\nLet's try to na\u00efvely make a brick-wall filter of length 0.1 sec, and look\nat the filter itself in the time domain and the frequency domain:\n\n"
      ], 
      "cell_type": "markdown", 
      "metadata": {}
    }, 
    {
      "execution_count": null, 
      "cell_type": "code", 
      "source": [
        "n = int(round(0.1 * sfreq)) + 1\nt = np.arange(-n // 2, n // 2) / sfreq  # center our sinc\nh = np.sinc(2 * f_p * t) / (4 * np.pi)\nplot_filter(h, sfreq, freq, gain, 'Sinc (0.1 sec)', flim=flim)"
      ], 
      "outputs": [], 
      "metadata": {
        "collapsed": false
      }
    }, 
    {
      "source": [
        "This is not so good! Making the filter 10 times longer (1 sec) gets us a\nbit better stop-band suppression, but still has a lot of ringing in\nthe time domain. Note the x-axis is an order of magnitude longer here,\nand the filter has a correspondingly much longer group delay (again equal\nto half the filter length, or 0.5 seconds):\n\n"
      ], 
      "cell_type": "markdown", 
      "metadata": {}
    }, 
    {
      "execution_count": null, 
      "cell_type": "code", 
      "source": [
        "n = int(round(1. * sfreq)) + 1\nt = np.arange(-n // 2, n // 2) / sfreq\nh = np.sinc(2 * f_p * t) / (4 * np.pi)\nplot_filter(h, sfreq, freq, gain, 'Sinc (1.0 sec)', flim=flim)"
      ], 
      "outputs": [], 
      "metadata": {
        "collapsed": false
      }
    }, 
    {
      "source": [
        "Let's make the stop-band tighter still with a longer filter (10 sec),\nwith a resulting larger x-axis:\n\n"
      ], 
      "cell_type": "markdown", 
      "metadata": {}
    }, 
    {
      "execution_count": null, 
      "cell_type": "code", 
      "source": [
        "n = int(round(10. * sfreq)) + 1\nt = np.arange(-n // 2, n // 2) / sfreq\nh = np.sinc(2 * f_p * t) / (4 * np.pi)\nplot_filter(h, sfreq, freq, gain, 'Sinc (10.0 sec)', flim=flim)"
      ], 
      "outputs": [], 
      "metadata": {
        "collapsed": false
      }
    }, 
    {
      "source": [
        "Now we have very sharp frequency suppression, but our filter rings for the\nentire second. So this na\u00efve method is probably not a good way to build\nour low-pass filter.\n\nFortunately, there are multiple established methods to design FIR filters\nbased on desired response characteristics. These include:\n\n    1. The Remez_ algorithm (:func:`scipy.signal.remez`, `MATLAB firpm`_)\n    2. Windowed FIR design (:func:`scipy.signal.firwin2`, `MATLAB fir2`_)\n    3. Least squares designs (:func:`scipy.signal.firls`, `MATLAB firls`_)\n    4. Frequency-domain design (construct filter in Fourier\n       domain and use an :func:`IFFT <scipy.fftpack.ifft>` to invert it)\n\n<div class=\"alert alert-info\"><h4>Note</h4><p>Remez and least squares designs have advantages when there are\n          \"do not care\" regions in our frequency response. However, we want\n          well controlled responses in all frequency regions.\n          Frequency-domain construction is good when an arbitrary response\n          is desired, but generally less clean (due to sampling issues) than\n          a windowed approach for more straightfroward filter applications.\n          Since our filters (low-pass, high-pass, band-pass, band-stop)\n          are fairly simple and we require precisel control of all frequency\n          regions, here we will use and explore primarily windowed FIR\n          design.</p></div>\n\nIf we relax our frequency-domain filter requirements a little bit, we can\nuse these functions to construct a lowpass filter that instead has a\n*transition band*, or a region between the pass frequency $f_p$\nand stop frequency $f_s$, e.g.:\n\n"
      ], 
      "cell_type": "markdown", 
      "metadata": {}
    }, 
    {
      "execution_count": null, 
      "cell_type": "code", 
      "source": [
        "trans_bandwidth = 10  # 10 Hz transition band\nf_s = f_p + trans_bandwidth  # = 50 Hz\n\nfreq = [0., f_p, f_s, nyq]\ngain = [1., 1., 0., 0.]\nax = plt.subplots(1, figsize=third_height)[1]\ntitle = '%s Hz lowpass with a %s Hz transition' % (f_p, trans_bandwidth)\nplot_ideal_filter(freq, gain, ax, title=title, flim=flim)"
      ], 
      "outputs": [], 
      "metadata": {
        "collapsed": false
      }
    }, 
    {
      "source": [
        "Accepting a shallower roll-off of the filter in the frequency domain makes\nour time-domain response potentially much better. We end up with a\nsmoother slope through the transition region, but a *much* cleaner time\ndomain signal. Here again for the 1 sec filter:\n\n"
      ], 
      "cell_type": "markdown", 
      "metadata": {}
    }, 
    {
      "execution_count": null, 
      "cell_type": "code", 
      "source": [
        "h = signal.firwin2(n, freq, gain, nyq=nyq)\nplot_filter(h, sfreq, freq, gain, 'Windowed 10-Hz transition (1.0 sec)',\n            flim=flim)"
      ], 
      "outputs": [], 
      "metadata": {
        "collapsed": false
      }
    }, 
    {
      "source": [
        "Since our lowpass is around 40 Hz with a 10 Hz transition, we can actually\nuse a shorter filter (5 cycles at 10 Hz = 0.5 sec) and still get okay\nstop-band attenuation:\n\n"
      ], 
      "cell_type": "markdown", 
      "metadata": {}
    }, 
    {
      "execution_count": null, 
      "cell_type": "code", 
      "source": [
        "n = int(round(sfreq * 0.5)) + 1\nh = signal.firwin2(n, freq, gain, nyq=nyq)\nplot_filter(h, sfreq, freq, gain, 'Windowed 10-Hz transition (0.5 sec)',\n            flim=flim)"
      ], 
      "outputs": [], 
      "metadata": {
        "collapsed": false
      }
    }, 
    {
      "source": [
        "But then if we shorten the filter too much (2 cycles of 10 Hz = 0.2 sec),\nour effective stop frequency gets pushed out past 60 Hz:\n\n"
      ], 
      "cell_type": "markdown", 
      "metadata": {}
    }, 
    {
      "execution_count": null, 
      "cell_type": "code", 
      "source": [
        "n = int(round(sfreq * 0.2)) + 1\nh = signal.firwin2(n, freq, gain, nyq=nyq)\nplot_filter(h, sfreq, freq, gain, 'Windowed 10-Hz transition (0.2 sec)',\n            flim=flim)"
      ], 
      "outputs": [], 
      "metadata": {
        "collapsed": false
      }
    }, 
    {
      "source": [
        "If we want a filter that is only 0.1 seconds long, we should probably use\nsomething more like a 25 Hz transition band (0.2 sec = 5 cycles @ 25 Hz):\n\n"
      ], 
      "cell_type": "markdown", 
      "metadata": {}
    }, 
    {
      "execution_count": null, 
      "cell_type": "code", 
      "source": [
        "trans_bandwidth = 25\nf_s = f_p + trans_bandwidth\nfreq = [0, f_p, f_s, nyq]\nh = signal.firwin2(n, freq, gain, nyq=nyq)\nplot_filter(h, sfreq, freq, gain, 'Windowed 50-Hz transition (0.2 sec)',\n            flim=flim)"
      ], 
      "outputs": [], 
      "metadata": {
        "collapsed": false
      }
    }, 
    {
      "source": [
        "So far we have only discussed *acausal* filtering, which means that each\nsample at each time point $t$ is filtered using samples that come\nafter ($t + \\Delta t$) *and* before ($t - \\Delta t$) $t$.\nIn this sense, each sample is influenced by samples that come both before\nand after it. This is useful in many cases, espcially because it does not\ndelay the timing of events.\n\nHowever, sometimes it can be beneficial to use *causal* filtering,\nwhereby each sample $t$ is filtered only using time points that came\nafter it.\n\nNote that the delay is variable (whereas for linear/zero-phase filters it\nis constant) but small in the pass-band. Unlike zero-phase filters, which\nrequire time-shifting backward the output of a linear-phase filtering stage\n(and thus becoming acausal), minimum-phase filters do not require any\ncompensation to achieve small delays in the passband. Note that as an\nartifact of the minimum phase filter construction step, the filter does\nnot end up being as steep as the linear/zero-phase version.\n\nWe can construct a minimum-phase filter from our existing linear-phase\nfilter with the ``minimum_phase`` function (that will be in SciPy 0.19's\n:mod:`scipy.signal`), and note that the falloff is not as steep:\n\n"
      ], 
      "cell_type": "markdown", 
      "metadata": {}
    }, 
    {
      "execution_count": null, 
      "cell_type": "code", 
      "source": [
        "h_min = mne.fixes.minimum_phase(h)\nplot_filter(h_min, sfreq, freq, gain, 'Minimum-phase', flim=flim)"
      ], 
      "outputs": [], 
      "metadata": {
        "collapsed": false
      }
    }, 
    {
      "source": [
        "\nApplying FIR filters\n--------------------\n\nNow lets look at some practical effects of these filters by applying\nthem to some data.\n\nLet's construct a Gaussian-windowed sinusoid (i.e., Morlet imaginary part)\nplus noise (random + line). Note that the original, clean signal contains\nfrequency content in both the pass band and transition bands of our\nlow-pass filter.\n\n"
      ], 
      "cell_type": "markdown", 
      "metadata": {}
    }, 
    {
      "execution_count": null, 
      "cell_type": "code", 
      "source": [
        "dur = 10.\ncenter = 2.\nmorlet_freq = f_p\ntlim = [center - 0.2, center + 0.2]\ntticks = [tlim[0], center, tlim[1]]\nflim = [20, 70]\n\nx = np.zeros(int(sfreq * dur) + 1)\nblip = morlet(sfreq, [morlet_freq], n_cycles=7)[0].imag / 20.\nn_onset = int(center * sfreq) - len(blip) // 2\nx[n_onset:n_onset + len(blip)] += blip\nx_orig = x.copy()\n\nrng = np.random.RandomState(0)\nx += rng.randn(len(x)) / 1000.\nx += np.sin(2. * np.pi * 60. * np.arange(len(x)) / sfreq) / 2000."
      ], 
      "outputs": [], 
      "metadata": {
        "collapsed": false
      }
    }, 
    {
      "source": [
        "Filter it with a shallow cutoff, linear-phase FIR (which allows us to\ncompensate for the constant filter delay):\n\n"
      ], 
      "cell_type": "markdown", 
      "metadata": {}
    }, 
    {
      "execution_count": null, 
      "cell_type": "code", 
      "source": [
        "transition_band = 0.25 * f_p\nf_s = f_p + transition_band\nfilter_dur = 6.6 / transition_band  # sec\nn = int(sfreq * filter_dur)\nfreq = [0., f_p, f_s, sfreq / 2.]\ngain = [1., 1., 0., 0.]\n# This would be equivalent:\n# h = signal.firwin2(n, freq, gain, nyq=sfreq / 2.)\nh = mne.filter.create_filter(x, sfreq, l_freq=None, h_freq=f_p)\nx_shallow = np.convolve(h, x)[len(h) // 2:]\n\nplot_filter(h, sfreq, freq, gain, 'MNE-Python 0.14 default', flim=flim)"
      ], 
      "outputs": [], 
      "metadata": {
        "collapsed": false
      }
    }, 
    {
      "source": [
        "This is actually set to become the default type of filter used in MNE-Python\nin 0.14 (see `tut_filtering_in_python`).\n\nLet's also filter with the MNE-Python 0.13 default, which is a\nlong-duration, steep cutoff FIR that gets applied twice:\n\n"
      ], 
      "cell_type": "markdown", 
      "metadata": {}
    }, 
    {
      "execution_count": null, 
      "cell_type": "code", 
      "source": [
        "transition_band = 0.5  # Hz\nf_s = f_p + transition_band\nfilter_dur = 10.  # sec\nn = int(sfreq * filter_dur)\nfreq = [0., f_p, f_s, sfreq / 2.]\ngain = [1., 1., 0., 0.]\n# This would be equivalent\n# h = signal.firwin2(n, freq, gain, nyq=sfreq / 2.)\nh = mne.filter.create_filter(x, sfreq, l_freq=None, h_freq=f_p,\n                             h_trans_bandwidth=transition_band,\n                             filter_length='%ss' % filter_dur)\nx_steep = np.convolve(np.convolve(h, x)[::-1], h)[::-1][len(h) - 1:-len(h) - 1]\n\nplot_filter(h, sfreq, freq, gain, 'MNE-Python 0.13 default', flim=flim)"
      ], 
      "outputs": [], 
      "metadata": {
        "collapsed": false
      }
    }, 
    {
      "source": [
        "Let's also filter it with the MNE-C default, which is a long-duration\nsteep-slope FIR filter designed using frequency-domain techniques:\n\n"
      ], 
      "cell_type": "markdown", 
      "metadata": {}
    }, 
    {
      "execution_count": null, 
      "cell_type": "code", 
      "source": [
        "h = mne.filter.design_mne_c_filter(sfreq, l_freq=None, h_freq=f_p + 2.5)\nx_mne_c = np.convolve(h, x)[len(h) // 2:]\n\ntransition_band = 5  # Hz (default in MNE-C)\nf_s = f_p + transition_band\nfreq = [0., f_p, f_s, sfreq / 2.]\ngain = [1., 1., 0., 0.]\nplot_filter(h, sfreq, freq, gain, 'MNE-C default', flim=flim)"
      ], 
      "outputs": [], 
      "metadata": {
        "collapsed": false
      }
    }, 
    {
      "source": [
        "And now an example of a minimum-phase filter:\n\n"
      ], 
      "cell_type": "markdown", 
      "metadata": {}
    }, 
    {
      "execution_count": null, 
      "cell_type": "code", 
      "source": [
        "h = mne.filter.create_filter(x, sfreq, l_freq=None, h_freq=f_p,\n                             phase='minimum')\nx_min = np.convolve(h, x)\ntransition_band = 0.25 * f_p\nf_s = f_p + transition_band\nfilter_dur = 6.6 / transition_band  # sec\nn = int(sfreq * filter_dur)\nfreq = [0., f_p, f_s, sfreq / 2.]\ngain = [1., 1., 0., 0.]\nplot_filter(h, sfreq, freq, gain, 'Minimum-phase filter', flim=flim)"
      ], 
      "outputs": [], 
      "metadata": {
        "collapsed": false
      }
    }, 
    {
      "source": [
        "Both the MNE-Python 0.13 and MNE-C filhters have excellent frequency\nattenuation, but it comes at a cost of potential\nringing (long-lasting ripples) in the time domain. Ringing can occur with\nsteep filters, especially on signals with frequency content around the\ntransition band. Our Morlet wavelet signal has power in our transition band,\nand the time-domain ringing is thus more pronounced for the steep-slope,\nlong-duration filter than the shorter, shallower-slope filter:\n\n"
      ], 
      "cell_type": "markdown", 
      "metadata": {}
    }, 
    {
      "execution_count": null, 
      "cell_type": "code", 
      "source": [
        "axes = plt.subplots(1, 2)[1]\n\n\ndef plot_signal(x, offset):\n    t = np.arange(len(x)) / sfreq\n    axes[0].plot(t, x + offset)\n    axes[0].set(xlabel='Time (sec)', xlim=t[[0, -1]])\n    X = fftpack.fft(x)\n    freqs = fftpack.fftfreq(len(x), 1. / sfreq)\n    mask = freqs >= 0\n    X = X[mask]\n    freqs = freqs[mask]\n    axes[1].plot(freqs, 20 * np.log10(np.abs(X)))\n    axes[1].set(xlim=flim)\n\nyticks = np.arange(6) / -30.\nyticklabels = ['Original', 'Noisy', 'FIR-shallow (0.14)', 'FIR-steep (0.13)',\n               'FIR-steep (MNE-C)', 'Minimum-phase']\nplot_signal(x_orig, offset=yticks[0])\nplot_signal(x, offset=yticks[1])\nplot_signal(x_shallow, offset=yticks[2])\nplot_signal(x_steep, offset=yticks[3])\nplot_signal(x_mne_c, offset=yticks[4])\nplot_signal(x_min, offset=yticks[5])\naxes[0].set(xlim=tlim, title='FIR, Lowpass=%d Hz' % f_p, xticks=tticks,\n            ylim=[-0.200, 0.025], yticks=yticks, yticklabels=yticklabels,)\nfor text in axes[0].get_yticklabels():\n    text.set(rotation=45, size=8)\naxes[1].set(xlim=flim, ylim=(-60, 10), xlabel='Frequency (Hz)',\n            ylabel='Magnitude (dB)')\nmne.viz.tight_layout()\nplt.show()"
      ], 
      "outputs": [], 
      "metadata": {
        "collapsed": false
      }
    }, 
    {
      "source": [
        "IIR filters\n===========\n\nMNE-Python also offers IIR filtering functionality that is based on the\nmethods from :mod:`scipy.signal`. Specifically, we use the general-purpose\nfunctions :func:`scipy.signal.iirfilter` and :func:`scipy.signal.iirdesign`,\nwhich provide unified interfaces to IIR filter design.\n\nDesigning IIR filters\n---------------------\n\nLet's continue with our design of a 40 Hz low-pass filter, and look at\nsome trade-offs of different IIR filters.\n\nOften the default IIR filter is a `Butterworth filter`_, which is designed\nto have a *maximally flat pass-band*. Let's look at a few orders of filter,\ni.e., a few different number of coefficients used and therefore steepness\nof the filter:\n\n<div class=\"alert alert-info\"><h4>Note</h4><p>Notice that the group delay (which is related to the phase) of\n          the IIR filters below are not constant. In the FIR case, we can\n          design so-called linear-phase filters that have a constant group\n          delay, and thus compensate for the delay (making the filter\n          acausal) if necessary. This cannot be done with IIR filters, as\n          they have a non-linear phase (non-constant group delay). As the\n          filter order increases, the phase distortion near and in the\n          transition band worsens. However, if acausal (forward-backward)\n          filtering can be used, e.g. with :func:`scipy.signal.filtfilt`,\n          these phase issues can theoretically be mitigated.</p></div>\n\n"
      ], 
      "cell_type": "markdown", 
      "metadata": {}
    }, 
    {
      "execution_count": null, 
      "cell_type": "code", 
      "source": [
        "sos = signal.iirfilter(2, f_p / nyq, btype='low', ftype='butter', output='sos')\nplot_filter(dict(sos=sos), sfreq, freq, gain, 'Butterworth order=2', flim=flim)\n\n# Eventually this will just be from scipy signal.sosfiltfilt, but 0.18 is\n# not widely adopted yet (as of June 2016), so we use our wrapper...\nsosfiltfilt = mne.fixes.get_sosfiltfilt()\nx_shallow = sosfiltfilt(sos, x)"
      ], 
      "outputs": [], 
      "metadata": {
        "collapsed": false
      }
    }, 
    {
      "source": [
        "The falloff of this filter is not very steep.\n\n<div class=\"alert alert-info\"><h4>Note</h4><p>Here we have made use of second-order sections (SOS)\n          by using :func:`scipy.signal.sosfilt` and, under the\n          hood, :func:`scipy.signal.zpk2sos` when passing the\n          ``output='sos'`` keyword argument to\n          :func:`scipy.signal.iirfilter`. The filter definitions\n          given in tut_filtering_basics_ use the polynomial\n          numerator/denominator (sometimes called \"tf\") form ``(b, a)``,\n          which are theoretically equivalent to the SOS form used here.\n          In practice, however, the SOS form can give much better results\n          due to issues with numerical precision (see\n          :func:`scipy.signal.sosfilt` for an example), so SOS should be\n          used when possible to do IIR filtering.</p></div>\n\nLet's increase the order, and note that now we have better attenuation,\nwith a longer impulse response:\n\n"
      ], 
      "cell_type": "markdown", 
      "metadata": {}
    }, 
    {
      "execution_count": null, 
      "cell_type": "code", 
      "source": [
        "sos = signal.iirfilter(8, f_p / nyq, btype='low', ftype='butter', output='sos')\nplot_filter(dict(sos=sos), sfreq, freq, gain, 'Butterworth order=8', flim=flim)\nx_steep = sosfiltfilt(sos, x)"
      ], 
      "outputs": [], 
      "metadata": {
        "collapsed": false
      }
    }, 
    {
      "source": [
        "There are other types of IIR filters that we can use. For a complete list,\ncheck out the documentation for :func:`scipy.signal.iirdesign`. Let's\ntry a Chebychev (type I) filter, which trades off ripple in the pass-band\nto get better attenuation in the stop-band:\n\n"
      ], 
      "cell_type": "markdown", 
      "metadata": {}
    }, 
    {
      "execution_count": null, 
      "cell_type": "code", 
      "source": [
        "sos = signal.iirfilter(8, f_p / nyq, btype='low', ftype='cheby1', output='sos',\n                       rp=1)  # dB of acceptable pass-band ripple\nplot_filter(dict(sos=sos), sfreq, freq, gain,\n            'Chebychev-1 order=8, ripple=1 dB', flim=flim)"
      ], 
      "outputs": [], 
      "metadata": {
        "collapsed": false
      }
    }, 
    {
      "source": [
        "And if we can live with even more ripple, we can get it slightly steeper,\nbut the impulse response begins to ring substantially longer (note the\ndifferent x-axis scale):\n\n"
      ], 
      "cell_type": "markdown", 
      "metadata": {}
    }, 
    {
      "execution_count": null, 
      "cell_type": "code", 
      "source": [
        "sos = signal.iirfilter(8, f_p / nyq, btype='low', ftype='cheby1', output='sos',\n                       rp=6)\nplot_filter(dict(sos=sos), sfreq, freq, gain,\n            'Chebychev-1 order=8, ripple=6 dB', flim=flim)"
      ], 
      "outputs": [], 
      "metadata": {
        "collapsed": false
      }
    }, 
    {
      "source": [
        "Applying IIR filters\n--------------------\n\nNow let's look at how our shallow and steep Butterworth IIR filters\nperform on our Morlet signal from before:\n\n"
      ], 
      "cell_type": "markdown", 
      "metadata": {}
    }, 
    {
      "execution_count": null, 
      "cell_type": "code", 
      "source": [
        "axes = plt.subplots(1, 2)[1]\nyticks = np.arange(4) / -30.\nyticklabels = ['Original', 'Noisy', 'Butterworth-2', 'Butterworth-8']\nplot_signal(x_orig, offset=yticks[0])\nplot_signal(x, offset=yticks[1])\nplot_signal(x_shallow, offset=yticks[2])\nplot_signal(x_steep, offset=yticks[3])\naxes[0].set(xlim=tlim, title='IIR, Lowpass=%d Hz' % f_p, xticks=tticks,\n            ylim=[-0.125, 0.025], yticks=yticks, yticklabels=yticklabels,)\nfor text in axes[0].get_yticklabels():\n    text.set(rotation=45, size=8)\naxes[1].set(xlim=flim, ylim=(-60, 10), xlabel='Frequency (Hz)',\n            ylabel='Magnitude (dB)')\nmne.viz.adjust_axes(axes)\nmne.viz.tight_layout()\nplt.show()"
      ], 
      "outputs": [], 
      "metadata": {
        "collapsed": false
      }
    }, 
    {
      "source": [
        "Some pitfalls of filtering\n==========================\n\nMultiple recent papers have noted potential risks of drawing\nerrant inferences due to misapplication of filters.\n\nLow-pass problems\n-----------------\n\nFilters in general, especially those that are acausal (zero-phase), can make\nactivity appear to occur earlier or later than it truly did. As\nmentioned in VanRullen 2011 [3]_, investigations of commonly (at the time)\nused low-pass filters created artifacts when they were applied to smulated\ndata. However, such deleterious effects were minimal in many real-world\nexamples in Rousselet 2012 [5]_.\n\nPerhaps more revealing, it was noted in Widmann & Schr\u00f6ger 2012 [6]_ that\nthe problematic low-pass filters from VanRullen 2011 [3]_:\n\n   1. Used a least-squares design (like :func:`scipy.signal.firls`) that\n      included \"do-not-care\" transition regions, which can lead to\n      uncontrolled behavior.\n   2. Had a filter length that was independent of the transition bandwidth,\n      which can cause excessive ringing and signal distortion.\n\n\nHigh-pass problems\n------------------\n\nWhen it comes to high-pass filtering, using corner frequencies above 0.1 Hz\nwere found in Acunzo *et al.* 2012 [4]_ to:\n\n   \"...generate a systematic bias easily leading to misinterpretations of\n   neural activity.\u201d\n\nIn a related paper, Widmann *et al.* 2015 [7]_ also came to suggest a 0.1 Hz\nhighpass. And more evidence followed in Tanner *et al.* 2015 [8]_ of such\ndistortions. Using data from language ERP studies of semantic and syntactic\nprocessing (i.e., N400 and P600), using a high-pass above 0.3 Hz caused\nsignificant effects to be introduced implausibly early when compared to the\nunfiltered data. From this, the authors suggested the optimal high-pass\nvalue for language processing to be 0.1 Hz.\n\nWe can recreate a problematic simulation from Tanner *et al.* 2015 [8]_:\n\n   \"The simulated component is a single-cycle cosine wave with an amplitude\n   of 5\u00b5V, onset of 500 ms poststimulus, and duration of 800 ms. The\n   simulated component was embedded in 20 s of zero values to avoid\n   filtering edge effects... Distortions [were] caused by 2 Hz low-pass and\n   high-pass filters... No visible distortion to the original waveform\n   [occurred] with 30 Hz low-pass and 0.01 Hz high-pass filters...\n   Filter frequencies correspond to the half-amplitude (-6 dB) cutoff\n   (12 dB/octave roll-off).\"\n\n<div class=\"alert alert-info\"><h4>Note</h4><p>This simulated signal contains energy not just within the\n          pass-band, but also within the transition and stop-bands -- perhaps\n          most easily understood because the signal has a non-zero DC value,\n          but also because it is a shifted cosine that has been\n          *windowed* (here multiplied by a rectangular window), which\n          makes the cosine and DC frequencies spread to other frequencies\n          (multiplication in time is convolution in frequency, so multiplying\n          by a rectangular window in the time domain means convolving a sinc\n          function with the impulses at DC and the cosine frequency in the\n          frequency domain).</p></div>\n\n\n"
      ], 
      "cell_type": "markdown", 
      "metadata": {}
    }, 
    {
      "execution_count": null, 
      "cell_type": "code", 
      "source": [
        "x = np.zeros(int(2 * sfreq))\nt = np.arange(0, len(x)) / sfreq - 0.2\nonset = np.where(t >= 0.5)[0][0]\ncos_t = np.arange(0, int(sfreq * 0.8)) / sfreq\nsig = 2.5 - 2.5 * np.cos(2 * np.pi * (1. / 0.8) * cos_t)\nx[onset:onset + len(sig)] = sig\n\niir_lp_30 = signal.iirfilter(2, 30. / sfreq, btype='lowpass')\niir_hp_p1 = signal.iirfilter(2, 0.1 / sfreq, btype='highpass')\niir_lp_2 = signal.iirfilter(2, 2. / sfreq, btype='lowpass')\niir_hp_2 = signal.iirfilter(2, 2. / sfreq, btype='highpass')\nx_lp_30 = signal.filtfilt(iir_lp_30[0], iir_lp_30[1], x, padlen=0)\nx_hp_p1 = signal.filtfilt(iir_hp_p1[0], iir_hp_p1[1], x, padlen=0)\nx_lp_2 = signal.filtfilt(iir_lp_2[0], iir_lp_2[1], x, padlen=0)\nx_hp_2 = signal.filtfilt(iir_hp_2[0], iir_hp_2[1], x, padlen=0)\n\nxlim = t[[0, -1]]\nylim = [-2, 6]\nxlabel = 'Time (sec)'\nylabel = 'Amplitude ($\\mu$V)'\ntticks = [0, 0.5, 1.3, t[-1]]\naxes = plt.subplots(2, 2)[1].ravel()\nfor ax, x_f, title in zip(axes, [x_lp_2, x_lp_30, x_hp_2, x_hp_p1],\n                          ['LP$_2$', 'LP$_{30}$', 'HP$_2$', 'LP$_{0.1}$']):\n    ax.plot(t, x, color='0.5')\n    ax.plot(t, x_f, color='k', linestyle='--')\n    ax.set(ylim=ylim, xlim=xlim, xticks=tticks,\n           title=title, xlabel=xlabel, ylabel=ylabel)\nmne.viz.adjust_axes(axes)\nmne.viz.tight_layout()\nplt.show()"
      ], 
      "outputs": [], 
      "metadata": {
        "collapsed": false
      }
    }, 
    {
      "source": [
        "Similarly, in a P300 paradigm reported by Kappenman & Luck 2010 [12]_,\nthey found that applying a 1 Hz high-pass decreased the probaility of\nfinding a significant difference in the N100 response, likely because\nthe P300 response was smeared (and inverted) in time by the high-pass\nfilter such that it tended to cancel out the increased N100. However,\nthey nonetheless note that some high-passing can still be useful to deal\nwith drifts in the data.\n\nEven though these papers generally advise a 0.1 HZ or lower frequency for\na high-pass, it is important to keep in mind (as most authors note) that\nfiltering choices should depend on the frequency content of both the\nsignal(s) of interest and the noise to be suppressed. For example, in\nsome of the MNE-Python examples involving `ch_sample_data`,\nhigh-pass values of around 1 Hz are used when looking at auditory\nor visual N100 responses, because we analyze standard (not deviant) trials\nand thus expect that contamination by later or slower components will\nbe limited.\n\nBaseline problems (or solutions?)\n---------------------------------\n\nIn an evolving discussion, Tanner *et al.* 2015 [8]_ suggest using baseline\ncorrection to remove slow drifts in data. However, Maess *et al.* 2016 [9]_\nsuggest that baseline correction, which is a form of high-passing, does\nnot offer substantial advantages over standard high-pass filtering.\nTanner *et al.* [10]_ rebutted that baseline correction can correct for\nproblems with filtering.\n\nTo see what they mean, consider again our old simulated signal ``x`` from\nbefore:\n\n"
      ], 
      "cell_type": "markdown", 
      "metadata": {}
    }, 
    {
      "execution_count": null, 
      "cell_type": "code", 
      "source": [
        "def baseline_plot(x):\n    all_axes = plt.subplots(3, 2)[1]\n    for ri, (axes, freq) in enumerate(zip(all_axes, [0.1, 0.3, 0.5])):\n        for ci, ax in enumerate(axes):\n            if ci == 0:\n                iir_hp = signal.iirfilter(4, freq / sfreq, btype='highpass',\n                                          output='sos')\n                x_hp = sosfiltfilt(iir_hp, x, padlen=0)\n            else:\n                x_hp -= x_hp[t < 0].mean()\n            ax.plot(t, x, color='0.5')\n            ax.plot(t, x_hp, color='k', linestyle='--')\n            if ri == 0:\n                ax.set(title=('No ' if ci == 0 else '') +\n                       'Baseline Correction')\n            ax.set(xticks=tticks, ylim=ylim, xlim=xlim, xlabel=xlabel)\n            ax.set_ylabel('%0.1f Hz' % freq, rotation=0,\n                          horizontalalignment='right')\n        mne.viz.adjust_axes(axes)\n    mne.viz.tight_layout()\n    plt.suptitle(title)\n    plt.show()\n\nbaseline_plot(x)"
      ], 
      "outputs": [], 
      "metadata": {
        "collapsed": false
      }
    }, 
    {
      "source": [
        "In respose, Maess *et al.* 2016 [11]_ note that these simulations do not\naddress cases of pre-stimulus activity that is shared across conditions, as\napplying baseline correction will effectively copy the topology outside the\nbaseline period. We can see this if we give our signal ``x`` with some\nconsistent pre-stimulus activity, which makes everything look bad.\n\n<div class=\"alert alert-info\"><h4>Note</h4><p>An important thing to keep in mind with these plots is that they\n          are for a single simulated sensor. In multielectrode recordings\n          the topology (i.e., spatial pattiern) of the pre-stimulus activity\n          will leak into the post-stimulus period. This will likely create a\n          spatially varying distortion of the time-domain signals, as the\n          averaged pre-stimulus spatial pattern gets subtracted from the\n          sensor time courses.</p></div>\n\nPutting some activity in the baseline period:\n\n"
      ], 
      "cell_type": "markdown", 
      "metadata": {}
    }, 
    {
      "execution_count": null, 
      "cell_type": "code", 
      "source": [
        "n_pre = (t < 0).sum()\nsig_pre = 1 - np.cos(2 * np.pi * np.arange(n_pre) / (0.5 * n_pre))\nx[:n_pre] += sig_pre\nbaseline_plot(x)"
      ], 
      "outputs": [], 
      "metadata": {
        "collapsed": false
      }
    }, 
    {
      "source": [
        "Both groups seem to acknowledge that the choices of filtering cutoffs, and\nperhaps even the application of baseline correction, depend on the\ncharacteristics of the data being investigated, especially when it comes to:\n\n   1. The frequency content of the underlying evoked activity relative\n      to the filtering parameters.\n   2. The validity of the assumption of no consistent evoked activity\n      in the baseline period.\n\nWe thus recommend carefully applying baseline correction and/or high-pass\nvalues based on the characteristics of the data to be analyzed.\n\n\nFiltering defaults\n==================\n\n\nDefaults in MNE-Python\n----------------------\n\nMost often, filtering in MNE-Python is done at the :class:`mne.io.Raw` level,\nand thus :func:`mne.io.Raw.filter` is used. This function under the hood\n(among other things) calls :func:`mne.filter.filter_data` to actually\nfilter the data, which by default applies a zero-phase FIR filter designed\nusing :func:`scipy.signal.firwin2`. In Widmann *et al.* 2015 [7]_, they\nsuggest a specific set of parameters to use for high-pass filtering,\nincluding:\n\n    \"... providing a transition bandwidth of 25% of the lower passband\n    edge but, where possible, not lower than 2 Hz and otherwise the\n    distance from the passband edge to the critical frequency.\u201d\n\nIn practice, this means that for each high-pass value ``l_freq`` or\nlow-pass value ``h_freq`` below, you would get this corresponding\n``l_trans_bandwidth`` or ``h_trans_bandwidth``, respectively,\nif the sample rate were 100 Hz (i.e., Nyquist frequency of 50 Hz):\n\n+------------------+-------------------+-------------------+\n| l_freq or h_freq | l_trans_bandwidth | h_trans_bandwidth |\n+==================+===================+===================+\n|             0.01 |              0.01 |               2.0 |\n+------------------+-------------------+-------------------+\n|              0.1 |               0.1 |               2.0 |\n+------------------+-------------------+-------------------+\n|              1.0 |               1.0 |               2.0 |\n+------------------+-------------------+-------------------+\n|              2.0 |               2.0 |               2.0 |\n+------------------+-------------------+-------------------+\n|              4.0 |               2.0 |               2.0 |\n+------------------+-------------------+-------------------+\n|              8.0 |               2.0 |               2.0 |\n+------------------+-------------------+-------------------+\n|             10.0 |               2.5 |               2.5 |\n+------------------+-------------------+-------------------+\n|             20.0 |               5.0 |               5.0 |\n+------------------+-------------------+-------------------+\n|             40.0 |              10.0 |              10.0 |\n+------------------+-------------------+-------------------+\n|             45.0 |             11.25 |              5.0  |\n+------------------+-------------------+-------------------+\n|             48.0 |              12.0 |              2.0  |\n+------------------+-------------------+-------------------+\n\nMNE-Python has adopted this definition for its high-pass (and low-pass)\ntransition bandwidth choices when using ``l_trans_bandwidth='auto'`` and\n``h_trans_bandwidth='auto'``.\n\nTo choose the filter length automatically with ``filter_length='auto'``,\nthe reciprocal of the shortest transition bandwidth is used to ensure\ndecent attenuation at the stop frequency. Specifically, the reciprocal\n(in samples) is multiplied by 6.2, 6.6, or 11.0 for the Hann, Hamming,\nor Blackman windows, respectively as selected by the ``fir_window``\nargument.\n\n<div class=\"alert alert-info\"><h4>Note</h4><p>These multiplicative factors are double what is given in\n          Ifeachor and Jervis [2]_ (p. 357). The window functions have a\n          smearing effect on the frequency response; I&J thus take the\n          approach of setting the stop frequency as\n          $f_s = f_p + f_{trans} / 2.$, but our stated definitions of\n          $f_s$ and $f_{trans}$ do not\n          allow us to do this in a nice way. Instead, we increase our filter\n          length to achieve acceptable (20+ dB) attenuation by\n          $f_s = f_p + f_{trans}$, and excellent (50+ dB)\n          attenuation by $f_s + f_{trans}$ (and usually earlier).</p></div>\n\nIn 0.14, we default to using a Hamming window in filter design, as it\nprovides up to 53 dB of stop-band attenuation with small pass-band ripple.\n\n<div class=\"alert alert-info\"><h4>Note</h4><p>In band-pass applications, often a low-pass filter can operate\n          effectively with fewer samples than the high-pass filter, so\n          it is advisable to apply the high-pass and low-pass separately.</p></div>\n\nFor more information on how to use the\nMNE-Python filtering functions with real data, consult the preprocessing\ntutorial on `tut_artifacts_filter`.\n\nDefaults in MNE-C\n-----------------\nMNE-C by default uses:\n\n   1. 5 Hz transition band for low-pass filters.\n   2. 3-sample transition band for high-pass filters.\n   3. Filter length of 8197 samples.\n\nThe filter is designed in the frequency domain, creating a linear-phase\nfilter such that the delay is compensated for as is done with the MNE-Python\n``phase='zero'`` filtering option.\n\nSquared-cosine ramps are used in the transition regions. Because these\nare used in place of more gradual (e.g., linear) transitions,\na given transition width will result in more temporal ringing but also more\nrapid attenuation than the same transition width in windowed FIR designs.\n\nThe default filter length will generally have excellent attenuation\nbut long ringing for the sample rates typically encountered in M-EEG data\n(e.g. 500-2000 Hz).\n\nDefaults in other software\n--------------------------\nA good but possibly outdated comparison of filtering in various software\npackages is available in [7]_. Briefly:\n\n* EEGLAB\n    MNE-Python in 0.14 defaults to behavior very similar to that of EEGLAB,\n    see the `EEGLAB filtering FAQ`_ for more information.\n* Fieldrip\n    By default FieldTrip applies a forward-backward Butterworth IIR filter\n    of order 4 (band-pass and band-stop filters) or 2 (for low-pass and\n    high-pass filters). Similar filters can be achieved in MNE-Python when\n    filtering with :meth:`raw.filter(..., method='iir') <mne.io.Raw.filter>`\n    (see also :func:`mne.filter.construct_iir_filter` for options).\n    For more inforamtion, see e.g. `FieldTrip band-pass documentation`_.\n\nSummary\n=======\n\nWhen filtering, there are always tradeoffs that should be considered.\nOne important tradeoff is between time-domain characteristics (like ringing)\nand frequency-domain attenuation characteristics (like effective transition\nbandwidth). Filters with sharp frequency cutoffs can produce outputs that\nring for a long time when they operate on signals with frequency content\nin the transition band. In general, therefore, the wider a transition band\nthat can be tolerated, the better behaved the filter will be in the time\ndomain.\n\nReferences\n==========\n\n.. [1] Parks TW, Burrus CS (1987). Digital Filter Design.\n       New York: Wiley-Interscience.\n.. [2] Ifeachor, E. C., & Jervis, B. W. (2002). Digital Signal Processing:\n       A Practical Approach. Prentice Hall.\n.. [3] Vanrullen, R. (2011). Four common conceptual fallacies in mapping\n       the time course of recognition. Perception Science, 2, 365.\n.. [4] Acunzo, D. J., MacKenzie, G., & van Rossum, M. C. W. (2012).\n       Systematic biases in early ERP and ERF components as a result\n       of high-pass filtering. Journal of Neuroscience Methods,\n       209(1), 212\u2013218. http://doi.org/10.1016/j.jneumeth.2012.06.011\n.. [5] Rousselet, G. A. (2012). Does filtering preclude us from studying\n       ERP time-courses? Frontiers in Psychology, 3(131)\n.. [6] Widmann, A., & Schr\u00f6ger, E. (2012). Filter effects and filter\n       artifacts in the analysis of electrophysiological data.\n       Perception Science, 233.\n.. [7] Widmann, A., Schr\u00f6ger, E., & Maess, B. (2015). Digital filter\n       design for electrophysiological data \u2013 a practical approach.\n       Journal of Neuroscience Methods, 250, 34\u201346.\n.. [8] Tanner, D., Morgan-Short, K., & Luck, S. J. (2015).\n       How inappropriate high-pass filters can produce artifactual effects\n       and incorrect conclusions in ERP studies of language and cognition.\n       Psychophysiology, 52(8), 997\u20131009. http://doi.org/10.1111/psyp.12437\n.. [9] Maess, B., Schr\u00f6ger, E., & Widmann, A. (2016).\n       High-pass filters and baseline correction in M/EEG analysis.\n       Commentary on: \u201cHow inappropriate high-pass filters can produce\n       artefacts and incorrect conclusions in ERP studies of language\n       and cognition.\u201d Journal of Neuroscience Methods, 266, 164\u2013165.\n.. [10] Tanner, D., Norton, J. J. S., Morgan-Short, K., & Luck, S. J. (2016).\n       On high-pass filter artifacts (they\u2019re real) and baseline correction\n       (it\u2019s a good idea) in ERP/ERMF analysis.\n.. [11] Maess, B., Schr\u00f6ger, E., & Widmann, A. (2016).\n       High-pass filters and baseline correction in M/EEG analysis-continued\n       discussion. Journal of Neuroscience Methods, 266, 171\u2013172.\n       Journal of Neuroscience Methods, 266, 166\u2013170.\n.. [12] Kappenman E. & Luck, S. (2010). The effects of impedance on data\n       quality and statistical significance in ERP recordings.\n       Psychophysiology, 47, 888-904.\n\n\n"
      ], 
      "cell_type": "markdown", 
      "metadata": {}
    }
  ], 
  "metadata": {
    "kernelspec": {
      "display_name": "Python 2", 
      "name": "python2", 
      "language": "python"
    }, 
    "language_info": {
      "mimetype": "text/x-python", 
      "nbconvert_exporter": "python", 
      "name": "python", 
      "file_extension": ".py", 
      "version": "2.7.13", 
      "pygments_lexer": "ipython2", 
      "codemirror_mode": {
        "version": 2, 
        "name": "ipython"
      }
    }
  }
}