<jupyter_start><jupyter_text>Herding Spikes 2 detection (with locally exclusive spikes)recall lower than expected? bandpass filter greatly improves it<jupyter_code>detect_peaks_threshold_list = [15, 20, 25]
radius_size = 50 #um
time_err_threshold = 0.0002 # 0.2 ms
precisions = []
recalls = []
accuracys = []

filtered_recording = st.bandpass_filter(recording, freq_min=300.0, freq_max=6000.0)
print(f"search radius size: {radius_size}")
for threshold in detect_peaks_threshold_list:
    with silence():
        peaks = scp.hs2_spike_detection(recording=filtered_recording, to_localize=False, detect_threshold=threshold)
    marked = {}   
    accuracy, precision, recall = evaluate_detection(recording=recording, detected_spikes=peaks, gt_sorting=sorting_gt, time_error_threshold=time_err_threshold, radius_size=radius_size)
    print(f"detection threshold: {threshold}, accuracy: {accuracy}, precision: {precision}, recall: {recall}\n")
    precisions.append(precision)
    recalls.append(recall)
    accuracys.append(accuracy)

plot_metrics(detect_peaks_threshold_list, accuracys, precisions, recalls)<jupyter_output>search radius size: 50
total_detected_spikes: 771, total_gt_spikes: 747, matched_spikes: 720
detection threshold: 15, accuracy: 0.9022556390977443, precision: 0.933852140077821, recall: 0.963855421686747

total_detected_spikes: 751, total_gt_spikes: 747, matched_spikes: 715
detection threshold: 20, accuracy: 0.913154533844189, precision: 0.952063914780293, recall: 0.9571619812583668

total_detected_spikes: 715, total_gt_spikes: 747, matched_spikes: 689
detection threshold: 25, accuracy: 0.8913324708926261, precision: 0.9636363636363636, recall: 0.9223560910307899<jupyter_text>detect_peaks (method = by channel)<jupyter_code>detect_peaks_threshold_list = [1, 3, 5, 7, 9] # In median absolute deviations (MAD).
radius_size = 0 #um
time_err_threshold = 0.0002 # 0.2 ms
precisions = []
recalls = []
accuracys = []

print(f"search radius size: {radius_size}, time_err_threshold: {time_err_threshold * 1000} ms")
for threshold in detect_peaks_threshold_list:
    peaks = scp.detect_peaks(recording=recording, detect_threshold=threshold)
    marked = {}   
    accuracy, precision, recall = evaluate_detection(recording=recording, detected_spikes=peaks, gt_sorting=sorting_gt, time_error_threshold=time_err_threshold, radius_size=radius_size)
    print(f"detection threshold: {threshold}, accuracy: {accuracy}, precision: {precision}, recall: {recall}\n")
    precisions.append(precision)
    recalls.append(recall)
    accuracys.append(accuracy)

plot_metrics(detect_peaks_threshold_list, accuracys, precisions, recalls)
detect_peaks_threshold_list = [1, 3, 5, 7, 9] # In median absolute deviations (MAD).
radius_size = 50 #um
time_err_threshold = 0.0002 # 0.2 ms
precisions = []
recalls = []
accuracys = []

print(f"search radius size: {radius_size}, time_err_threshold: {time_err_threshold * 1000} ms")
for threshold in detect_peaks_threshold_list:
    peaks = scp.detect_peaks(recording=recording, detect_threshold=threshold)
    marked = {}   
    accuracy, precision, recall = evaluate_detection(recording=recording, detected_spikes=peaks, gt_sorting=sorting_gt, time_error_threshold=time_err_threshold, radius_size=radius_size)
    print(f"detection threshold: {threshold}, accuracy: {accuracy}, precision: {precision}, recall: {recall}\n")
    precisions.append(precision)
    recalls.append(recall)
    accuracys.append(accuracy)

plot_metrics(detect_peaks_threshold_list, accuracys, precisions, recalls)<jupyter_output>search radius size: 50, time_err_threshold: 0.2 ms
total_detected_spikes: 851133, total_gt_spikes: 747, matched_spikes: 747
detection threshold: 1, accuracy: 0.0008776536686980766, precision: 0.0008776536686980766, recall: 1.0

total_detected_spikes: 24501, total_gt_spikes: 747, matched_spikes: 747
detection threshold: 3, accuracy: 0.03048855148769438, precision: 0.03048855148769438, recall: 1.0

total_detected_spikes: 7054, total_gt_spikes: 747, matched_spikes: 747
detection threshold: 5, accuracy: 0.10589736319818542, precision: 0.10589736319818542, recall: 1.0

total_detected_spikes: 3453, total_gt_spikes: 747, matched_spikes: 728
detection threshold: 7, accuracy: 0.20967741935483872, precision: 0.21083116130900667, recall: 0.9745649263721553

total_detected_spikes: 1646, total_gt_spikes: 747, matched_spikes: 423
detection threshold: 9, accuracy: 0.2147208121827411, precision: 0.25698663426488455, recall: 0.5662650602409639<jupyter_text>detect_peaks (method = locally exclusive)<jupyter_code>detect_peaks_threshold_list = [1, 3, 5, 7, 9] # In median absolute deviations (MAD).
radius_size = 50 #um
time_err_threshold = 0.0002 # 0.2 ms
precisions = []
recalls = []
accuracys = []

print(f"search radius size: {radius_size}, time_err_threshold: {time_err_threshold * 1000} ms")
for threshold in detect_peaks_threshold_list:
    peaks = scp.detect_peaks(recording=recording, detect_threshold=threshold, method='locally_exclusive')
    marked = {}   
    accuracy, precision, recall = evaluate_detection(recording=recording, detected_spikes=peaks, gt_sorting=sorting_gt, time_error_threshold=time_err_threshold, radius_size=radius_size)
    print(f"detection threshold: {threshold}, accuracy: {accuracy}, precision: {precision}, recall: {recall}\n")
    precisions.append(precision)
    recalls.append(recall)
    accuracys.append(accuracy)

plot_metrics(detect_peaks_threshold_list, accuracys, precisions, recalls)<jupyter_output>search radius size: 50, time_err_threshold: 0.2 ms
total_detected_spikes: 204279, total_gt_spikes: 747, matched_spikes: 747
detection threshold: 1, accuracy: 0.00365676354397662, precision: 0.00365676354397662, recall: 1.0

total_detected_spikes: 9465, total_gt_spikes: 747, matched_spikes: 747
detection threshold: 3, accuracy: 0.07892234548335975, precision: 0.07892234548335975, recall: 1.0

total_detected_spikes: 1628, total_gt_spikes: 747, matched_spikes: 747
detection threshold: 5, accuracy: 0.45884520884520885, precision: 0.45884520884520885, recall: 1.0

total_detected_spikes: 1050, total_gt_spikes: 747, matched_spikes: 728
detection threshold: 7, accuracy: 0.6810102899906455, precision: 0.6933333333333334, recall: 0.9745649263721553

total_detected_spikes: 446, total_gt_spikes: 747, matched_spikes: 423
detection threshold: 9, accuracy: 0.5493506493506494, precision: 0.9484304932735426, recall: 0.5662650602409639<jupyter_text>detect_peaks evaluation by making spikes locally exclusive<jupyter_code>def radius_search_2(radius_size, detected_spikes, peak_channel_ind, spike_to_match, time_error_threshold):
    duplicate_spikes = []
    channel_inds_in_radius = get_radius_channels(radius_size, channel_ind=peak_channel_ind)
    channel_inds_in_radius.append(peak_channel_ind)
    # print(f"channel inds in radius: {channel_inds_in_radius}")

    for channel_ind in channel_inds_in_radius:
        channel_spikes = detected_spikes[detected_spikes['channel_ind'] == channel_ind]
        # channel_spikes = np.asarray(sorted(channel_spikes, key=lambda x: x['sample_ind']))
        if len(channel_spikes) == 0: continue
        channel_spikes = channel_spikes['sample_ind'] / sampling_freq
        for spike_ind, spike_time in enumerate(channel_spikes):
            if (channel_ind, spike_ind) in marked: 
                # print("already marked")
                continue
            if abs(spike_time - spike_to_match) <= time_error_threshold:
                duplicate_spikes.append((channel_ind, spike_ind))
    
    if len(duplicate_spikes) == 0: return [(None, None)]
    else: return duplicate_spikes
def evaluate_detection_2(recording, detected_spikes, gt_sorting, time_error_threshold, radius_size):
    matched_gt_spikes = 0
    total_gt_spikes = 0
    for unit in sorting_gt.get_unit_ids():
        peak_channel_ind = int(extremum_channels_inds[unit])
        gt_spikes = sorting_gt.get_unit_spike_train(unit) / sampling_freq
        total_gt_spikes += len(gt_spikes)

        # print(f"unit: {unit}, num gt spikes: {len(gt_spikes)}")
        for gt_ind, spike_time in enumerate(gt_spikes):
            channel_spike_ind_list = radius_search_2(radius_size=radius_size, detected_spikes=peaks, peak_channel_ind=peak_channel_ind, spike_to_match=spike_time, time_error_threshold=time_error_threshold)
            # print(f"channel_spike_ind_list: {channel_spike_ind_list}")
            # print(f"finished looking for gt spike with ind: {gt_ind}")

            matched_gt_spike = False
            for (channel_ind, spike_ind) in channel_spike_ind_list:
                if channel_ind != None and spike_ind != None:
                    marked[(channel_ind, spike_ind)] = True
                    matched_gt_spike = True
            if matched_gt_spike == True: matched_gt_spikes += 1
    recall = matched_gt_spikes / (matched_gt_spikes + (total_gt_spikes - matched_gt_spikes))
    precision = matched_gt_spikes / (matched_gt_spikes + (len(detected_spikes) - len(marked)))
    accuracy = matched_gt_spikes / (matched_gt_spikes + (total_gt_spikes - matched_gt_spikes) + (len(detected_spikes) - len(marked)))
    print(f"total_detected_spikes: {len(detected_spikes)}, matched_detected_spikes: {len(marked)}, total_gt_spikes: {total_gt_spikes}, matched_gt_spikes: {matched_gt_spikes}")
    return accuracy, precision, recall
detect_peaks_threshold_list = [3, 5, 7, 9] # In median absolute deviations (MAD).
radius_size = 50 #um
time_err_threshold = 0.0004
precisions = []
recalls = []
accuracys = []

print(f"search radius size: {radius_size}, time_err_threshold: {time_err_threshold * 1000} ms")
for threshold in detect_peaks_threshold_list:
    peaks = scp.detect_peaks(recording=filtered_recording, detect_threshold=threshold)
    marked = {}   
    accuracy, precision, recall = evaluate_detection_2(recording=recording, detected_spikes=peaks, gt_sorting=sorting_gt, time_error_threshold=time_err_threshold, radius_size=radius_size)
    print(f"detection threshold: {threshold}, accuracy: {accuracy}, precision: {precision}, recall: {recall}\n")
    precisions.append(precision)
    recalls.append(recall)
    accuracys.append(accuracy)

plot_metrics(detect_peaks_threshold_list, accuracys, precisions, recalls)<jupyter_output>search radius size: 50, time_err_threshold: 0.4 ms
total_detected_spikes: 18195, matched_detected_spikes: 7679, total_gt_spikes: 747, matched_gt_spikes: 747
detection threshold: 3, accuracy: 0.06632335967326644, precision: 0.06632335967326644, recall: 1.0

total_detected_spikes: 4918, matched_detected_spikes: 4563, total_gt_spikes: 747, matched_gt_spikes: 746
detection threshold: 5, accuracy: 0.676950998185118, precision: 0.6775658492279746, recall: 0.998661311914324

total_detected_spikes: 2673, matched_detected_spikes: 2648, total_gt_spikes: 747, matched_gt_spikes: 742
detection threshold: 7, accuracy: 0.961139896373057, precision: 0.9674054758800521, recall: 0.9933065595716198

total_detected_spikes: 1521, matched_detected_spikes: 1521, total_gt_spikes: 747, matched_gt_spikes: 526
detection threshold: 9, accuracy: 0.7041499330655957, precision: 1.0, recall: 0.7041499330655957<jupyter_text>Visualising accuracy per unit and SNR<jupyter_code>def evaluate_detection_per_unit(recording, detected_spikes, gt_sorting, time_error_threshold, radius_size):
    precision_per_unit = []
    recalls_per_unit = []
    accuracy_per_unit = []

    total_peak_channel_detected_spikes = 0
    for unit in sorting_gt.get_unit_ids():
        matched_gt_spikes = 0
        marked = {}
        peak_channel_ind = int(extremum_channels_inds[unit])
        gt_spikes = sorting_gt.get_unit_spike_train(unit) / sampling_freq
        crt_unit_gt_spikes = len(gt_spikes)

        for spike_time in gt_spikes:
            channel_ind, spike_ind = radius_search(radius_size=radius_size, detected_spikes=peaks, peak_channel_ind=peak_channel_ind, spike_to_match=spike_time, time_error_threshold=time_error_threshold)
            if channel_ind != None and spike_ind != None:
                marked[(channel_ind, spike_ind)] = True
                matched_gt_spikes += 1

        peak_channel_spikes = detected_spikes[detected_spikes['channel_ind'] == peak_channel_ind]
        total_peak_channel_detected_spikes += len(peak_channel_spikes)
        print(f"spikes detected on peak channel: {len(peak_channel_spikes)}, matched gt spikes: {matched_gt_spikes}, crt gt unit spikes: {crt_unit_gt_spikes}")
        recall = matched_gt_spikes / (matched_gt_spikes + (crt_unit_gt_spikes - matched_gt_spikes))
        precision = matched_gt_spikes / (matched_gt_spikes + (len(peak_channel_spikes) - len(marked)))
        accuracy = matched_gt_spikes / (matched_gt_spikes + (crt_unit_gt_spikes - matched_gt_spikes) + (len(peak_channel_spikes) - len(marked)))

        precision_per_unit.append(precision)
        recalls_per_unit.append(recall)
        accuracy_per_unit.append(accuracy)

        print(f"unit ind: {unit}, SNR: {SNRs[unit]:.2f}, accuracy: {accuracy:.2f}, precision: {precision:.2f}, recall: {recall:.2f}\n")

    # print(f"total_detected_spikes: {len(detected_spikes)}, total_gt_spikes: {total_gt_spikes}, matched_spikes: {matched_gt_spikes}")
    # return matched_gt_spikes, total_gt_spikes, recall
    print(f"total spikes detected on peak channels associated with gt units: {total_peak_channel_detected_spikes}")
    print(f"total detected spikes amongst all channels: {len(detected_spikes)}")
    print(f"num of unique channels spikes were detected on: {len(np.unique(detected_spikes['channel_ind']))}")
    print(f"num of channels on the recording: {recording.get_num_channels()}")
    return accuracy_per_unit, precision_per_unit, recalls_per_unit
SNRs = st.compute_snrs(waveform_extractor=we)
print(SNRs)
radius_size = 0 #um
time_err_threshold = 0.0001
detect_threshold = 5

print(f"search radius size: {radius_size}, time_err_threshold: {time_err_threshold * 1000} ms")
peaks = scp.detect_peaks(recording=recording, detect_threshold=detect_threshold)
marked = {}   
accuracy_per_unit, precision_per_unit, recalls_per_unit = evaluate_detection_per_unit(recording=recording, detected_spikes=peaks, gt_sorting=sorting_gt, time_error_threshold=time_err_threshold, radius_size=radius_size)
# plot_metrics(list(SNRs.values()), accuracy_per_unit, precision_per_unit, recalls_per_unit, x_label='SNR')<jupyter_output>search radius size: 0, time_err_threshold: 0.1 ms
spikes detected on peak channel: 99, matched gt spikes: 53, crt gt unit spikes: 53
unit ind: #0, SNR: 23.45, accuracy: 0.54, precision: 0.54, recall: 1.00

spikes detected on peak channel: 129, matched gt spikes: 50, crt gt unit spikes: 50
unit ind: #1, SNR: 25.57, accuracy: 0.39, precision: 0.39, recall: 1.00

spikes detected on peak channel: 140, matched gt spikes: 42, crt gt unit spikes: 43
unit ind: #2, SNR: 12.00, accuracy: 0.30, precision: 0.30, recall: 0.98

spikes detected on peak channel: 234, matched gt spikes: 30, crt gt unit spikes: 30
unit ind: #3, SNR: 21.85, accuracy: 0.13, precision: 0.13, recall: 1.00

spikes detected on peak channel: 324, matched gt spikes: 47, crt gt unit spikes: 48
unit ind: #4, SNR: 6.90, accuracy: 0.14, precision: 0.15, recall: 0.98

spikes detected on peak channel: 297, matched gt spikes: 20, crt gt unit spikes: 37
unit ind: #5, SNR: 5.18, accuracy: 0.06, precision: 0.07, recall: 0.54

spikes dete[...]
{'fraction_non_alphanumeric': 0.09672376734188758, 'fraction_numerical': 0.0801146355761089, 'mean_word_length': 5.86058981233244, 'pattern_counts': {'":': 0, '<': 17, '<?xml version=': 0, '>': 16, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 1, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '5', 'dir': 'jupyter-structured-clean-dedup', 'id': '432792', 'license': 'MIT', 'n_tokens_mistral': 6173, 'n_tokens_neox': 5227, 'n_words': 1119, 'path': 'projects/detection-benchmark/detection_performance.ipynb', 'repo_name': 'oliche/spike-sorting-hackathon', 'size': '15353'}