Fitting scaling factors in waveform/template plots for irregular probe layouts

[zzhmark]

delta_y will be inversely weighted by the x distances to get a proper interval over y dimension. In this case contacts that are close vertically but far laterally will be less represented.

The issue:
#3745 (comment)

[alejoe91]

@zzhmark can you provide some screenshots before/after? :)

[zzhmark]

I've made further improvement for the code. This one uses Fourier Space to get the proper x and y interval. Here are the results on my probe.

before:

after:

The reason why the waveform scale was so small was because there was a tiny shift over very distant contacts in the same row/column.

I've also tested on generated probes like:

It's compatible with regular designs.

[zzhmark]

I removed a feature for thresholding the peaks in the power spectra. It was too high.

[alejoe91]

@zzhmark thanks for the implementation. Isn't it an overkill though? I think it will make it hard to maintain in the long term...

I'd rather go for an easier solution as you suggested here: 71821d0

[zzhmark]

Unfortunately, the earlier one didn't work well on one of my probes. The unique part is vulnerable to small numerical differences. I actually had an idea with medium complexity, which uses the distance between contacts as the weight, but it's not good for highly stretched layouts where deltax is far from deltay.

[alejoe91]

@zm711 @samuelgarcia what do you think?

I believe the implementation is too complex for a plotting function at this point, so I would rather go for a middle-way solution that gives OK results, but is not too complex

[zm711]

I would have to agree. I think we need to think about maintainability and this is super complicated. Absolutely agree that we need to improve this, but I think incremental fix would be better and as problem come up we get more complicated. But jumping to this level of complexity is tricky.

[zzhmark]

okay, I will get you a simpler one lately.

[zzhmark]

Let's look at this one, only 7 lines. To ensure the compatibility I set the minimum delta x and y as 20.

[zm711]

That implementation looks much better to me. I'm in the middle of some analyses, but I'll have to test this myself with a couple probe layouts and make sure it works.

[chrishalcrow]

This looks great @zzhmark ! I remember this going wrong for me a long time ago, and didn't know why!

Here's a before/after on a simple generated example:

import numpy as np
import spikeinterface.full as si
from probeinterface import generate_tetrode

tetrode = generate_tetrode()
tetrode.device_channel_indices = np.array([0,1,2,3])
rec, sort = si.generate_ground_truth_recording(probe = tetrode, num_channels=4)

sa = si.create_sorting_analyzer(recording=rec, sorting=sort)
sa.compute(["random_spikes", "noise_levels", "waveforms", "templates"])

BEFORE:

AFTER:

However, I think the defaults don't look good with Neuropixles. Here is before and after:

BEFORE:

AFTER

I suggest keeping the NP defaults, since this is a very common use case. And I think the tetrode data looks nice with them too. This is my suggestion:

delta_x = max(gap * scx / base, 32)
delta_y = max(gap * scy / base, 15)

which gives:

[zm711]

@chrishalcrow, don't you think the x on the tetrodes is just a hair too close? Again I don't think we will ever be perfect and I agree NP is more common so we want to make sure that always looks nice. I want to test one other probe (It has a v-shape so it's X is a bit weird so I'd like to test to see how this handles some of setups.

As I see it the normal geometries are
diamonds (eg. tetrodes)
squares/rectangle (eg. NP)
checkerboard (repeating diamonds really) (eg other NP)

then you could have some "weirder"

v -shaped
homemade

So I would love to just have an idea of how as many of these would look. I assume when Alessio/Sam made this originally it was with a specific subset of probes in mind (maybe some NPs)

[zzhmark]

@chrishalcrow @zm711

Unfortunately, I found this implementation is useless as in most cases min will be invoked. Now I've got you a slightly more complexed one that works well. It uses gaussian distribution to add more weights on nearby channels, and uses orientation between channels to cancel column neighbours' contributions to row interval and vice versa.

Plus, I'm not sure if it's appropriate to hard code an interval for NP, but considering the vertical gap of NP is only 10, I decreased the minimum in this version. Anyway, the automatic method can reliably yield a nice estimation whether it's checkerboard or grid or if any noise exists. In this way, the min value is most likely just for single row or column layout.

Some examples:

recording, sorting = generate_ground_truth_recording(
            durations=[30.0],
            sampling_frequency=28000.0,
            num_channels=32,
            num_units=10,
            generate_probe_kwargs=dict(
                num_columns=4,
                xpitch=20,
                ypitch=20,
                contact_shapes="circle",
                contact_shape_params={"radius": 6},
            ),
            generate_sorting_kwargs=dict(firing_rates=10.0, refractory_period_ms=4.0),
            noise_kwargs=dict(noise_levels=5.0, strategy="on_the_fly"),
            seed=2205,
        )

recording, sorting = generate_ground_truth_recording(
            durations=[30.0],
            sampling_frequency=28000.0,
            num_channels=32,
            num_units=10,
            generate_probe_kwargs=dict(
                num_columns=2,
                xpitch=20,
                ypitch=20,
                contact_shapes="circle",
                contact_shape_params={"radius": 6},
            ),
            generate_sorting_kwargs=dict(firing_rates=10.0, refractory_period_ms=4.0),
            noise_kwargs=dict(noise_levels=5.0, strategy="on_the_fly"),
            seed=2205,
        )

def generate_checkerboard(rows, cols, spacing=2, noise_level=1):
    """
    Generates 2D coordinates for a checkerboard layout with noise.
    
    :param rows: Number of rows
    :param cols: Number of columns
    :param spacing: Distance between points
    :param noise_level: Maximum noise to add (randomly between -noise_level and +noise_level)
    :return: List of (x, y) coordinates
    """
    coords = []
    
    for i in range(rows):
        for j in range(cols):
            # Checkerboard pattern offset: every other row shifts by half spacing
            x = j * spacing + (i % 2) * (spacing / 2)
            y = i * spacing
            
            # Add random noise within [-noise_level, +noise_level]
            x += np.random.uniform(-noise_level, noise_level)
            y += np.random.uniform(-noise_level, noise_level)
            
            coords.append((x, y))
    
    return coords

checkerboard_coords = generate_checkerboard(rows=8, cols=8, spacing=30, noise_level=1)

recording, sorting = generate_ground_truth_recording(
            durations=[30.0],
            sampling_frequency=28000.0,
            num_channels=64,
            num_units=10,
            generate_probe_kwargs=dict(
                num_columns=1,
                xpitch=20,
                ypitch=20,
                contact_shapes="circle",
                contact_shape_params={"radius": 6},
            ),
            generate_sorting_kwargs=dict(firing_rates=10.0, refractory_period_ms=4.0),
            noise_kwargs=dict(noise_levels=5.0, strategy="on_the_fly"),
            seed=2205,
        )

prb = recording.get_probe()
prb._contact_positions = np.array(checkerboard_coords)
recording = recording.set_probe(prb)