[feature request] - add causal backward filtering to correct hardware induced phase shift

[JuanPimientoCaicedo]

Hi, Guys.

I hope you are having a good time at the hackathon. I tried to find options to deal with the phase shift induced by the hardware filter in the Neuropixel probes. These are single pole RC filters, which can distort waveforms. I could not find any simple function to do that and I think the easiest way is to create a new filter class in the filter.py file.

I took the plunge and created it, and added some extra steps in the get_traces function to perform a backward filtering with either the scipy.signal.sosfilt or the scipy.signal.lfilter on the flipped array. I really don't know if that is enough to implement this new feature. But if you want to take a look at the changes, I can create a pull request.

Thank you!

[JoeZiminski]

Hey @JuanPimientoCaicedo cheers for this! others are better placed to comment but I can't off the top of my head think of any reason not to provide a causal filter as well as the zero-phase filters currently provided. I think it is no problem to reopen #2942 as a draft PR, it's useful as a place for discussion and also so people can give advice as development progresses.

Could you go into more detail on the phase-shift induced by the hardware filter in the NP probes? This may be useful for others to think about, I have not seen compensating for any filter-induced phase shift before in preprocessing steps for NP probes and so this might be a nice addition! Currently there is a phase_shift preprocessing function but that compensates for multiplexing time shift across channels ('phase_shift' here actually refers to how the time shift is implemented in frequency domain).

Some quick thoughts on the PR:

• At present the timeseries is flipped, I guess to induce the phase-shift in the direction such that it compensates for the phase-shift of the hardware filter of NP probes. To generalise this, it would be worth making this flip optional so the phase distortion direction can be set by the user.
• for checking against a bool, you can do if causal:. Python has some strange intricacies on equality checks, discussed here (that whole thread is good, but particularly that comment).

[JuanPimientoCaicedo]

Hi, @JoeZiminski, Thank you for your comments.

Right now, I am working with NHP Neuropixel probes 1.0, which divide the recording in two streams, the AP stream (with a high pass filter of 300 Hz) and and LF stream (with a bandpass filter of 0.5 to 500 Hz if I am not wrong). You can actually turn off the filter in the AP band, which then converts it in a broadband recording, but given that the amplitude range of the 1.0 version is +/- 1.2 mV, it has happened to us that ERPs saturate the signal and we experience clippings.

Hardware filters are causal filters by definition, meaning that the output is produced by past and present signals. This behavior induces phase shifts, which change across different frequencies. One way to deal with that is simply use another causal filter but in the flipped signal (that is what the scipy.signal.filtfilt does, forward-backward filtering). here you have a good explanation with examples: https://www.mathworks.com/help/signal/ref/filtfilt.html.

Thank you for the PR comments.

1. It makes sense, I don't know what could the reason be for someone to try to use a causal filter in the forward direction, I will reopen the PR and wait for other comments.
2. I will change the if causal == True: for if causal:.

[JoeZiminski]

Cool thanks @JuanPimientoCaicedo, I never thought about correcting for these phase shifts its a really cool idea. For reference, in the article linked in this post suggests in (at least rodent) NP1 probes at the AP band is bandpass-filtered with highpass 300 Hz then lowpass 10,000 Hz while the LF band is filtered lowpass at 1000 Hz. But I'm having a bit of trouble from that paper what is similar across NHP vs. rodent probes. Do you have a reference for those numbers? It would be great to collect as much info on this as possible, the changes in architecture across versions is also relevant for the multiplexing issue.

In the paper they say the measurable frequency range is 0.5 - 10,000 kHz, but then that the AP band is 0.3-10,000 kHz. I don't know enough about electronics to interpret the difference in the lower bound there. But I guess for the phase-shift the RC filters are the main thing. Do you know if the phase-response of a single-pole RC filter is fairly standard or whether it will vary a lot depending on the specific electronic components they have used? I guess the phase response from the single-pole RC filter can be fairly well emulated with a first-order digital (IIR in this case) filter?

[alejoe91]

HI guys, interesting topic! The phase-shift in this case is called "group delay" and it's the derivative of the phase response of the filter. For a single pole filter (like an RC filter), this should be a "bump" in correspondence of the cutoff frequency.

We can make some calculations on how much this delay would affect the signal..

Happy to hear what @oliche and @samuelgarcia think!

[oliche]

The full description of the filters is here: https://pubmed.ncbi.nlm.nih.gov/28422663/.

Those seem a bit more complicated than simple RC filters so I think the best is to ask manufacturer if they have the impulse response of those filters, as I'm a bit out of my depth with those buffered first order switched capacitor filters !

[JuanPimientoCaicedo]

Hi, thank you for your comments (I also saw your comment in the Neuropixels slack channel).

To continue adding to the discussion I would like to share this paper that explains how the filters affect the phase of the waveform, that is one of the reasons that I am interested about the filters: https://journals.plos.org/plosone/article?id=10.1371/journal.pone.0174790. If you are interested in the waveform shape (and be able to compare units recorded in Neuropixels with the filter ON or OFF), what I would do is to apply the reverse filter with coefficients produced by the digital Butterworth filter design:

That is what is already done in the Spikeinterface code,

             filter_coeff = scipy.signal.iirfilter(
                filter_order, band, fs=fs, analog=False, btype=btype, ftype=ftype, output=filter_mode
            )

Now, instead of using scipy.signal.sosfiltfilt, you use scipy.signal.sosfilt in the flipped numpy.flip signal (and then you have to flip it again).

causal filters have a nonlinear phase response, and my understanding is that they are well modelled by a butterworth filter,

https://www.sciencedirect.com/science/article/abs/pii/B9780124058668000140?via%3Dihub

I don't come from from an engineering background, so I am more than happy to have insight from you guys.

Now that you bring up the LFP, It really depends on your analysis questions. For example, if you want to know if a neuron is locked in phase with an underlying oscillation (for example, you are studying theta phase precession) It really does not matter whether your signal in the LFP stream has a phase shift or not (as long as you are aware you cannot make any claim about the peak or troughs or any specific section of the oscillation). Now, if you are analyzing the peaks and troughs themselves, or doing some phase amplitude coupling analysis. then it is important to backward filter those signals as well.

That means you should apply to your lf band a backward filter with the coefficients for a 1 kHz lowpass filter and apply two filters to your ap band:

1. a filter with a corner at 300 hz.
2. a filter with a corner at 1000 hz.

All the filters I am talking about are first order Butterworth filters.

Please, let me know whether I am wrong or not.

Cheers.

[jiumao2]

Hi @JuanPimientoCaicedo ,

I'm curious about the extent to which the phase shift might influence the results. While I've found that the impact on spike timing is generally negligible, I'm uncertain about its effect on spike sorting quality. Do you have any insights or data on how phase shift could affect spike sorting performance?

[JuanPimientoCaicedo]

Well, @jiumao2 , I don't think it would hugely affect spike sorting quality. Let us imagine that you are recording from a neuron that has an specific waveform shape, given that the filter is applied to the whole recording, it is probably safe to assume that all the individual waveforms of that unit are going to be distorted in the same way. If you generalize that logic, all the waveforms of the different units in the recording are going to be distorted linearly following the phase response of whatever filter you are using.

So, if I am right, the spike sorting should still be able to discriminate between waveforms between different units even when distorted.

I think the most important aspect of the phase shift caused by causal filtering is the fact that we do care about the waveform shape and we use it to classify cells. And in order to be able to compare different recording conditions (for example broadband recordings with highpass filtered data) you need to account for these phase differences.

here you have a reference that develops this more: https://journals.plos.org/plosone/article?id=10.1371/journal.pone.0174790

I am not totally sure, but I believe that if you apply this step before processing or in postprocessing, the difference should be negligible.

[zm711]

Code has been merged. Let's close this. Please open new issues if questions arise related to causal filtering.