<SYSTEM_TASK:>
Given the following text description, write Python code to implement the functionality described below step by step
<END_TASK>
<USER_TASK:>
Description:
 tridesclous example with locust dataset
Here a detail notebook that detail the locust dataset recodring by Christophe Pouzat.
This dataset is our classic.
It has be analyse yet by several tools in R, Python or C
Step1:  Download a small dataset
trideclous provide some datasets than can be downloaded with download_dataset.
Note this dataset contains 2 trials in 2 different files. (the original contains more!)
Each file is considers as a segment. tridesclous automatically deal with it.
Step2:  DataIO = define datasource and working dir
Theses 2 files are in RawData format this means binary format with interleaved channels.
Our dataset contains 2 segment of 28.8 second each, 4 channels. The sample rate is 15kHz.
Note that there is only one channel_group here (0).
Step3:  CatalogueConstructor
Step4:  Set some parameters
For a complet description of each params see main documentation.
Step5:  Estimate the median and mad of noiseon a small chunk of filtered signals.
This compute medians and mad of each channel.
Step6:  Run the main loop
Step7:  Clean peaks
Whis try to detect "bad peaks". They are artifact with very big amplitude value.
This peaks have to removed early and not be include in waveform extaction and pca.
Strange peak are tag with -9 (alien)
Step8:  sample some peaks for waveforms extraction
Take some waveforms in the signals n_left/n_right must be choosen carfully.
It is not necessary to intensive to select all peaks.
There are several method to select peaks the most simple is to select randomly.
Note that waveform are extracted now. It is too intensive. They are extacted on-the-fly when needed.
Step9:  Extact some noise snippet.
Here a step to extact snippet of noise (in between real peak)
Step10:  Project to smaller space
To reduce dimension of the waveforms (n_peaks, peak_width, n_channel) we chosse global_pca method which is appropriate for tetrode.
It consists of flatenning some_waveforms.shape (n_peaks, peak_width, n_channel) to (n_peaks, peak_width*n_channel) and then apply a standard PCA on it with sklearn.
Let's keep 5 component of it.
In case of more channel we could also do a 'by_channel_pca'
Step11:  find clusters
There are many option to cluster this features. here a simple one the well known kmeans method.
Unfortunatly we need to choose the number of cluster. Too bad... Let's take 12.
Later on we will be able to refine this manually.
Step12:  Manual clean with CatalogueWindow (or visual check)
This open a CatalogueWindow, here we can check, split merge, trash, play as long as we are not happy.
If we are happy, we can save the catalogue.
Step13:  Here a snappshot of CatalogueWindow
<img src="../doc/img/snapshot_cataloguewindow.png">
Auto clean of catatalogue
tridesclous offer some method for auto merge/trash/split some cluster.
After this we can re order cluster and construct the catalogue for the peeler.
Step14:  Peeler
Create and run the Peeler.
It should be pretty fast, here the computation take 1.32s for 28.8x2s of signal. This is a speed up of 43 over real time.
Step15:  Open PeelerWindow for visual checking
<ASSISTANT_TASK:>
Python Code:
%matplotlib inline
import time
import numpy as np
import matplotlib.pyplot as plt
import tridesclous as tdc
from tridesclous import DataIO, CatalogueConstructor, Peeler
Explanation: tridesclous example with locust dataset
Here a detail notebook that detail the locust dataset recodring by Christophe Pouzat.
This dataset is our classic.
It has be analyse yet by several tools in R, Python or C:
  * https://github.com/christophe-pouzat/PouzatDetorakisEuroScipy2014
  * https://github.com/christophe-pouzat/SortingABigDataSetWithPython
  * http://xtof.perso.math.cnrs.fr/locust.html
So we can compare the result.
The original datasets is here https://zenodo.org/record/21589
But we will work on a very small subset on github https://github.com/tridesclous/tridesclous_datasets/tree/master/locust
Overview
In tridesclous, the spike sorting is done in several step:
  * Define the datasource and working path. (class DataIO)
  * Construct a catalogue (class CatalogueConstructor) on a short chunk of data (for instance 60s)
    with several sub step :
    * signal pre-processing:
      * high pass filter (optional)
      * removal of common reference (optional)
      * noise estimation (median/mad) on a small chunk
      * normalisation = robust z-score
    * peak detection
    * select a subset of peaks. Unecessary and impossible to extract them all.
    * extract some waveform.
    * project theses waveforms in smaller dimention (pca, ...)
    * find cluster
    * auto clean cluster with euritics merge/split/trash
    * clean manually with GUI (class CatalogueWindow) : merge/split/trash
    * save centroids (median+mad + first and second derivative)
  * Apply the Peeler (class Peeler) on the long term signals. With several sub steps:
     * same signal preprocessing than before
     * find peaks
     * find the best cluster in catalogue for each peak
     * find the intersample jitter
     * remove the oversampled waveforms from the signals until there are not peaks in the signals.
     * check with GUI (class PeelerWindow)
End of explanation
#download dataset
localdir, filenames, params = tdc.download_dataset(name='locust')
print(filenames)
print(params)
Explanation: Download a small dataset
trideclous provide some datasets than can be downloaded with download_dataset.
Note this dataset contains 2 trials in 2 different files. (the original contains more!)
Each file is considers as a segment. tridesclous automatically deal with it.
End of explanation
#create a DataIO
import os, shutil
dirname = 'tridesclous_locust'
if os.path.exists(dirname):
    #remove is already exists
    shutil.rmtree(dirname)    
dataio = DataIO(dirname=dirname)
# feed DataIO
dataio.set_data_source(type='RawData', filenames=filenames, **params)
print(dataio)
#no need to setup the prb with dataio.set_probe_file() or dataio.download_probe()
#because it is a tetrode
    
Explanation: DataIO = define datasource and working dir
Theses 2 files are in RawData format this means binary format with interleaved channels.
Our dataset contains 2 segment of 28.8 second each, 4 channels. The sample rate is 15kHz.
Note that there is only one channel_group here (0).
End of explanation
cc = CatalogueConstructor(dataio=dataio)
print(cc)
Explanation: CatalogueConstructor
End of explanation
# global params
cc.set_global_params(chunksize=1024,mode='dense')
# pre processing filetring normalisation
cc.set_preprocessor_params(
            common_ref_removal=False,
            highpass_freq=300.,
            lowpass_freq=5000.,                                             
            lostfront_chunksize=64)
cc.set_peak_detector_params(
            peak_sign='-',
            relative_threshold=6.5,
            peak_span_ms=0.1)
Explanation: Set some parameters
For a complet description of each params see main documentation.
End of explanation
cc.estimate_signals_noise(seg_num=0, duration=15.)
print(cc.signals_medians)
print(cc.signals_mads)
Explanation: Estimate the median and mad of noiseon a small chunk of filtered signals.
This compute medians and mad of each channel.
End of explanation
t1 = time.perf_counter()
cc.run_signalprocessor(duration=60.)
t2 = time.perf_counter()
print('run_signalprocessor', t2-t1, 's')
print(cc)
Explanation: Run the main loop: signal preprocessing + peak detection
End of explanation
cc.clean_peaks(alien_value_threshold=60., mode='extremum_amplitude')
print(cc)
Explanation: Clean peaks
Whis try to detect "bad peaks". They are artifact with very big amplitude value.
This peaks have to removed early and not be include in waveform extaction and pca.
Strange peak are tag with -9 (alien)
End of explanation
cc.set_waveform_extractor_params(wf_left_ms=-1.5, wf_right_ms=2.5)
cc.sample_some_peaks(mode='rand', nb_max=20000)
Explanation: sample some peaks for waveforms extraction
Take some waveforms in the signals n_left/n_right must be choosen carfully.
It is not necessary to intensive to select all peaks.
There are several method to select peaks the most simple is to select randomly.
Note that waveform are extracted now. It is too intensive. They are extacted on-the-fly when needed.
End of explanation
cc.extract_some_noise(nb_snippet=300)
Explanation: Extact some noise snippet.
Here a step to extact snippet of noise (in between real peak)
End of explanation
cc.extract_some_features(method='global_pca', n_components=5)
print(cc)
Explanation: Project to smaller space
To reduce dimension of the waveforms (n_peaks, peak_width, n_channel) we chosse global_pca method which is appropriate for tetrode.
It consists of flatenning some_waveforms.shape (n_peaks, peak_width, n_channel) to (n_peaks, peak_width*n_channel) and then apply a standard PCA on it with sklearn.
Let's keep 5 component of it.
In case of more channel we could also do a 'by_channel_pca'
End of explanation
cc.find_clusters(method='kmeans', n_clusters=12)
print(cc)
Explanation: find clusters
There are many option to cluster this features. here a simple one the well known kmeans method.
Unfortunatly we need to choose the number of cluster. Too bad... Let's take 12.
Later on we will be able to refine this manually.
End of explanation
%gui qt5
import pyqtgraph as pg
app = pg.mkQApp()
win = tdc.CatalogueWindow(catalogueconstructor)
win.show()
app.exec_()    
Explanation: Manual clean with CatalogueWindow (or visual check)
This open a CatalogueWindow, here we can check, split merge, trash, play as long as we are not happy.
If we are happy, we can save the catalogue.
End of explanation
cc.auto_split_cluster()
    
cc.auto_merge_cluster()
    
cc.trash_low_extremum(min_extremum_amplitude=6.6)
cc.trash_small_cluster(minimum_size=10)
#order cluster by waveforms rms
cc.order_clusters(by='waveforms_rms')
#save the catalogue
cc.make_catalogue_for_peeler(inter_sample_oversampling=True)
Explanation: Here a snappshot of CatalogueWindow
<img src="../doc/img/snapshot_cataloguewindow.png">
Auto clean of catatalogue
tridesclous offer some method for auto merge/trash/split some cluster.
After this we can re order cluster and construct the catalogue for the peeler.
End of explanation
catalogue = dataio.load_catalogue(chan_grp=0)
peeler = Peeler(dataio)
peeler.change_params(catalogue=catalogue)
t1 = time.perf_counter()
peeler.run()
t2 = time.perf_counter()
print('peeler.run', t2-t1)
print()
for seg_num in range(dataio.nb_segment):
    spikes = dataio.get_spikes(seg_num)
    print('seg_num', seg_num, 'nb_spikes', spikes.size)
    
Explanation: Peeler
Create and run the Peeler.
It should be pretty fast, here the computation take 1.32s for 28.8x2s of signal. This is a speed up of 43 over real time.
End of explanation
%gui qt5
import pyqtgraph as pg
app = pg.mkQApp()
win = tdc.PeelerWindow(dataio=dataio, catalogue=initial_catalogue)
win.show()
app.exec_()
Explanation: Open PeelerWindow for visual checking
End of explanation
<END_TASK>