<jupyter_start><jupyter_text>Signal Quality EvaluationThis tutorial demonstrates how signal quality can be evaluated using MNE-NIRS... sidebar:: Relevant literature   <NAME> al., “PHOEBE: a method for real time mapping of   optodes-scalp coupling in functional near-infrared spectroscopy” in   Biomed. Opt. Express 7, 5104-5119 (2016).   <NAME>, Samuel, and <NAME>. "NIRSplot: a tool for   quality assessment of fNIRS scans." Optics and the Brain.   Optical Society of America, 2020.   <NAME>., <NAME>., <NAME>., <NAME>., <NAME>., &   <NAME>. (2014). Auditory cortex activation to natural speech and   simulated cochlear implant speech measured with functional near-infrared   spectroscopy. Hearing research, 309, 84-93.Ensuring your data is of high quality is essential to good scientific research.Evaluating data quality is an essential part of both data collection andanalysis.For assessing data quality during data acquisition tools such as PHOEBE arehighly recommended, these allow the researcher to address potential problemsimmediately and modify their experiment to improve the quality of recordeddata.It is also important to assess data quality during data analysis.MNE-Python and MNE-NIRS provides several mechanisms to allow researchersto evaluate the quality of their data and to include this information in theirdownstream processing.Dedicated tools exist for quality evaluatiuon such as <NAME> (2020).This tutorial demonstrates methods in MNE-NIRS and MNE-Python for determiningchannels with poor signalquality, and methods for determining time segments of data that are of lowquality in a subset of channels.Two methods are introduced in this tutorial. The scalp coupling index (SCI)and peak power (PP) metrics. Both of these methods examine the signal forthe presence of a heart beat signal, which indicates that the sensorswere in contact with the scalp. For further details see the papers listedin the relevant literature sidebar.<jupyter_code># sphinx_gallery_thumbnail_number = 7

# Authors: <NAME> <<EMAIL>>
#
# License: BSD (3-clause)

import os
import mne
import numpy as np
from itertools import compress
import matplotlib.pyplot as plt

from mne.preprocessing.nirs import optical_density
from mne_nirs.preprocessing import peak_power, scalp_coupling_index_windowed
from mne_nirs.visualisation import plot_timechannel_quality_metric<jupyter_output><empty_output><jupyter_text>Import dataHere we will work with the `fNIRS motor data `.We resample the data to make indexing exact times more convenient.We then convert the data to optical density and plot the raw signal.<jupyter_code>fnirs_data_folder = mne.datasets.fnirs_motor.data_path()
fnirs_cw_amplitude_dir = os.path.join(fnirs_data_folder, 'Participant-1')
raw_intensity = mne.io.read_raw_nirx(fnirs_cw_amplitude_dir, verbose=True)
raw_intensity.load_data().resample(4.0, npad="auto")
raw_od = optical_density(raw_intensity)
raw_od.plot(n_channels=55, duration=4000, show_scrollbars=False, clipping=None)<jupyter_output><empty_output><jupyter_text>From the above plot we observe that the data is relatively clean.Later we will add some artificial bad sections to demonstrate the qualityevaluation metrics.  Scalp Coupling IndexThe scalp coupling index (SCI) from Pollonini (2016) provides a measure ofthe quality of the signal for a channel over a specifiedmeasurement duration. See Pollonini (2016) for further details of thetheory and implementation. SCI evaluated over entire signalHere we calculate the SCI for each channel overthe entire signal and view the distribution of values.<jupyter_code>sci = mne.preprocessing.nirs.scalp_coupling_index(raw_od)
fig, ax = plt.subplots()
ax.hist(sci)
ax.set(xlabel='Scalp Coupling Index', ylabel='Count', xlim=[0, 1])<jupyter_output><empty_output><jupyter_text>We observe that most of the channels have a good SCI of 1, but a few channelshave a poorer score. We can list the channels with an SCI below a threshold.And we can mark these as bad in the MNE data. This way the functionsdownstream will know that the data in these channels is bad.For this example we set the threshold for a bad channel to be SCI < 0.8.We then print a list of the bad channels and observe their are 10 channels(five source-detector pairs) that are marked as bad.<jupyter_code>raw_od.info['bads'] = list(compress(raw_od.ch_names, sci < 0.7))
print(raw_od.info['bads'])<jupyter_output><empty_output><jupyter_text>We can plot the time course of the signal again and note that the badchannels are now displayed in grey to indicate they are bad.<jupyter_code>raw_od.plot(n_channels=55, duration=4000, show_scrollbars=False, clipping=None)<jupyter_output><empty_output><jupyter_text>Similarly, we can view the montage diagram and view where on the headthe bad channels were positioned.<jupyter_code>raw_od.plot_sensors()<jupyter_output><empty_output><jupyter_text>SCI evaluated over initial segment of signalThe scalp coupling index can be calculated over a limited section ofthe signal by cropping to the desired section. For example, if you wish toevaluate the data quality of the first 10 seconds of the signal.Note that the difference to evaluation of the entire signal was quite subtle,but this may vary depending on your experimental design and setup.<jupyter_code>sci = mne.preprocessing.nirs.scalp_coupling_index(raw_od.copy().crop(10))
fig, ax = plt.subplots()
ax.hist(sci)
ax.set(xlabel='Scalp Coupling Index', ylabel='Count', xlim=[0, 1])<jupyter_output><empty_output><jupyter_text>SCI evaluated over moving windowThe scalp coupling index can be calculated overwindowed chunks of the signal.This plot is based on the meg bad channel detection figuresavailable in mne-bids-pipeline.Black horizontal lines indicate channels that have been marked as bad(see above).The color in the left facet shows the raw scores,The color in the right facet indicates segments that are below the threshold.This is useful for determining if a channel becomes bad throughout anexperiment. This may occur due to movement dislodging the optode or manyother causes.In this example we set the time window to be 60 seconds, but the user maydefine a window length that is appropriate for the experiment.<jupyter_code>_, scores, times = scalp_coupling_index_windowed(raw_od, time_window=60)
plot_timechannel_quality_metric(raw_od, scores, times, threshold=0.7,
                                title="Scalp Coupling Index "
                                      "Quality Evaluation")<jupyter_output><empty_output><jupyter_text>Peak PowerIt may also be informative to view the quality of the signal at a finertime resolution. The Peak Power metric provides a quality metric evalautedover a 10 second window. This allows the user to view instances wherea subset of channels may be contaminated by artifacts for a short durationof the recording.<jupyter_code>raw_od, scores, times = peak_power(raw_od, time_window=10)
plot_timechannel_quality_metric(raw_od, scores, times, threshold=0.1,
                                title="Peak Power Quality Evaluation")<jupyter_output><empty_output><jupyter_text>Introduced NoiseThe above data is quite clean, so here we add some noise to channels todemonstrate that the algorithm is able to detect bad segments of data.We add a step like signal at two time instances to mimic an artifact.We add one artifact to channel S2-D4 at 1000 secondsand one artifact to channel S5-D13 at 2000 seconds.<jupyter_code># Add an artifact to channel S2-D4 at time 1000 seconds
raw_od._data[12, 4000:4080] = np.linspace(0, -0.5, 80) + raw_od._data[12, 4000]

# Add an artifact to channel S5-D13 at time 2000 seconds
raw_od._data[34, 8000:8080] = np.linspace(0, 0.5, 80) + raw_od._data[34, 8000]<jupyter_output><empty_output><jupyter_text>Next we plot just these channels to demonstrate that indeed an artifacthas been added.<jupyter_code>raw_od.copy().pick(picks=[12, 13, 34, 35]).\
    plot(n_channels=55, duration=40000, show_scrollbars=False,
         clipping=None, scalings={'fnirs_od': 0.2})<jupyter_output><empty_output><jupyter_text>Peak Power MetricTo determine specific time and channel instances where data is of lowquality the peak power metric can be utilise (Pollonini, 2016).Below the peak power is calculated and displayed for all channels.The default threshold of 0.1 is applied to the quality scores to determinebad segments of data.The left facet illustrates the raw quality scores.The right facet illustrates the channel and time segments that do not meetthe threshold criteria and are marked as bad.<jupyter_code>raw_od, scores, times = peak_power(raw_od, time_window=10)
plot_timechannel_quality_metric(raw_od, scores, times, threshold=0.1)<jupyter_output><empty_output><jupyter_text>AnnotationsSimilar to how entire channels were marked as bad above, the peak powerfunction annotates the raw data structure to indicate where the badsegments of data are. If we replot the subset of channels as above we notethat the bad segments are marked in red now, indicating this time sectioncontained bad data.Note that although the red line appears as if the bad marking is presentfor all channels, this is due to the plotting code. Internally the softwareknows that only certain channels had a bad segment, and downstream processingwill only treat the specified bad channels as bad. I.e. when generatingepochs for time 1000 s, no epochs would be generated for channel S2-D4, butchannels would be generated for S5-D13 (as the artifact was only presenton S2-D4).<jupyter_code>raw_od.copy().pick(picks=[12, 13, 34, 35]).\
    plot(n_channels=55, duration=40000, show_scrollbars=False,
         clipping=None, scalings={'fnirs_od': 0.2})<jupyter_output><empty_output><jupyter_text>These channel and time specific annotations are used by downstreamprocessing. For example, when extracting epochs if a specific channelhas an annotation indicating a bad segment in a specific time window, thenthe epochs will be discarded for that specific channel.Finally, we display the entire time signal after both channel and timespecific quality evaluation. Note that bad channels are now greyed outand bad time segments are marked in red.<jupyter_code>raw_od.plot(n_channels=55, duration=4000, show_scrollbars=False, clipping=None)<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.06151896265964707, 'fraction_numerical': 0.020766305937408597, 'mean_word_length': 5.571428571428571, 'pattern_counts': {'":': 0, '<': 72, '<?xml version=': 0, '>': 70, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 1, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '14', 'dir': 'jupyter-structured-clean-dedup', 'id': '94527', 'license': 'BSD-3-Clause', 'n_tokens_mistral': 3010, 'n_tokens_neox': 2753, 'n_words': 1220, 'path': 'examples/mne-nirs-website-examples/auto_examples/general/plot_22_quality.ipynb', 'repo_name': 'rob-luke/mne-nirs-docker', 'size': '10342'}