How to use the spikeinterface.extractors function in spikeinterface

=======================

Spike sorters generally output a set of units with corresponding spike trains. The :code:`toolkit.postprocessing`
submodule allows to combine the :code:`RecordingExtractor` and the sorted :code:`SortingExtractor` objects to perform
further postprocessing.
"""

import matplotlib.pylab as plt

import spikeinterface.extractors as se
import spikeinterface.toolkit as st

##############################################################################
# First, let's create a toy example:

recording, sorting = se.example_datasets.toy_example(num_channels=4, duration=10, seed=0)

##############################################################################
# Assuming the :code:`sorting` is the output of a spike sorter, the
# :code:`postprocessing` module allows to extract all relevant information
# from the paired recording-sorting.


##############################################################################
# Compute spike waveforms
# --------------------------
# 
# Waveforms are extracted with the :code:`get_unit_waveforms` function by
# extracting snippets of the recordings when spikes are detected. When
# waveforms are extracted, the can be loaded in the :code:`SortingExtractor`
# object as features. The ms before and after the spike event can be
# chosen. Waveforms are returned as a list of np.arrays (n\_spikes,

'''
Comparison Widgets Gallery
===================================

Here is a gallery of all the available widgets using SortingExtractor objects.
'''

import spikeinterface.extractors as se
import spikeinterface.widgets as sw

##############################################################################
# First, let's create a toy example with the `extractors` module:

recording, sorting_true = se.example_datasets.toy_example(duration=10, num_channels=4, seed=0)

##############################################################################
# Let's run some spike sorting:

import spikeinterface.sorters as ss

sorting_MS4 = ss.run_mountainsort4(recording)
sorting_KL = ss.run_klusta(recording)

##############################################################################
# Widgets using SortingComparison
# ---------------------------------
#
# We can compare the spike sorting output to the ground-truth sorting :code:`sorting_true` using the
# :code:`comparison` module. :code:`comp_MS4` and :code:`comp_KL` are :code:`SortingComparison` objects

##############################################################################
# Installed Extractors
# --------------------
#
# To check which extractors are useable in a given python environment, one can print the installed recording extractor
# list and the installed sorting extractor list. An example from a newly installed miniconda3 environment is shown
# below.
#
# First, import the spikeextractors package:

import spikeinterface.extractors as se

##############################################################################
# Then you can check the installed RecordingExtractor list:

print(se.installed_recording_extractor_list)

##############################################################################
# and the installed SortingExtractors list:

print(se.installed_sorting_extractor_list)

##############################################################################
# When trying to use an extractor that has not been installed in your environment, an installation message will appear
# indicating which python packages must be installed as a prerequisite to using the extractor:

exdir_file = 'path_to_exdir_file'
try:
    recording = se.ExdirRecordingExtractor(exdir_file)
except Exception as e:
    print(e)

"""
Validation Tutorial
======================

After spike sorting, you might want to validate the goodness of the sorted units. This can be done using the
:code:`toolkit.validation` submodule, which computes several quality metrics of the sorted units.

"""

import spikeinterface.extractors as se
import spikeinterface.toolkit as st

##############################################################################
# First, let's create a toy example:

recording, sorting = se.example_datasets.toy_example(num_channels=4, duration=10, seed=0)

##############################################################################
# The :code:`toolkit.validation` submodule has a :code:`MetricCalculator` class that enables to compute metrics in a
# compact and easy way. You first need to instantiate a :code:`MetricCalculator` object with the
# :code:`SortingExtractor` and :code:`RecordingExtractor` objects.

mc = st.validation.MetricCalculator(sorting, recording)

##############################################################################
# You can then compute metrics as follows:

mc.compute_metrics()

##############################################################################
# This is the list of the computed metrics:
print(list(mc.get_metrics_dict().keys()))

import matplotlib.pyplot as plt
import seaborn as sns

import spikeinterface.extractors as se
import spikeinterface.widgets as sw
from spikeinterface.comparison import GroundTruthStudy

##############################################################################
#  Setup study folder and run all sorters
# ------------------------------------------------------
# 
#  We first generate the folder.
#  this can take some time because recordings are copied inside the folder.


rec0, gt_sorting0 = se.example_datasets.toy_example(num_channels=4, duration=10, seed=10)
rec1, gt_sorting1 = se.example_datasets.toy_example(num_channels=4, duration=10, seed=0)
gt_dict = {
    'rec0': (rec0, gt_sorting0),
    'rec1': (rec1, gt_sorting1),
}
study_folder = 'a_study_folder'
study = GroundTruthStudy.create(study_folder, gt_dict)

##############################################################################
#  Then just run all sorters on all recordings in one functions.

#  sorter_list = st.sorters.available_sorters() # this get all sorters.
sorter_list = ['klusta', 'tridesclous', 'mountainsort4']
study.run_sorters(sorter_list, mode="keep")

##############################################################################

sampling_frequency = 30000
duration = 20
num_timepoints = int(sampling_frequency * duration)
num_units = 4
num_events = 1000

##############################################################################
# We generate some random events.

times = np.int_(np.sort(np.random.uniform(0, num_timepoints, num_events)))
labels = np.random.randint(1, num_units + 1, size=num_events)

##############################################################################
# And instantiate a :code:`NumpyRecordingExtractor`:

sorting = se.NumpySortingExtractor()
sorting.set_times_labels(times=times, labels=labels)
sorting.set_sampling_frequency(sampling_frequency=sampling_frequency)

##############################################################################
# We can now print properties that the :code:`SortingExtractor` retrieves from the underlying sorted dataset.

print('Unit ids = {}'.format(sorting.get_unit_ids()))
st = sorting.get_unit_spike_train(unit_id=1)
print('Num. events for unit 1 = {}'.format(len(st)))
st1 = sorting.get_unit_spike_train(unit_id=1, start_frame=0, end_frame=30000)
print('Num. events for first second of unit 1 = {}'.format(len(st1)))

##############################################################################
# Some extractors also implement a :code:`write` function. We can for example save our newly created sorting into
# MDA format (Mountainsort4 format):

sorting.set_sampling_frequency(sampling_frequency=sampling_frequency)

##############################################################################
# We can now print properties that the :code:`SortingExtractor` retrieves from the underlying sorted dataset.

print('Unit ids = {}'.format(sorting.get_unit_ids()))
st = sorting.get_unit_spike_train(unit_id=1)
print('Num. events for unit 1 = {}'.format(len(st)))
st1 = sorting.get_unit_spike_train(unit_id=1, start_frame=0, end_frame=30000)
print('Num. events for first second of unit 1 = {}'.format(len(st1)))

##############################################################################
# Some extractors also implement a :code:`write` function. We can for example save our newly created sorting into
# MDA format (Mountainsort4 format):

se.MdaSortingExtractor.write_sorting(sorting=sorting, save_path='firings_true.mda')

##############################################################################
# and read it back with the proper extractor:

sorting2 = se.MdaSortingExtractor(file_path='firings_true.mda',
                                  sampling_frequency=sampling_frequency)
print('Unit ids = {}'.format(sorting2.get_unit_ids()))
st = sorting2.get_unit_spike_train(unit_id=1)
print('Num. events for unit 1 = {}'.format(len(st)))
st1 = sorting2.get_unit_spike_train(unit_id=1, start_frame=0, end_frame=30000)
print('Num. events for first second of unit 1 = {}'.format(len(st1)))

##############################################################################
# Unit properties are name value pairs that we can store for any unit. We will now calculate a unit
# property and store it in the :code:`SortingExtractor`