Question about new bandpass filter implementation

[Edouard2laire]

Hello there,

As I was preparing the NIRSTORM course, I found something weird with the new bandpass filter implementation (2019). In our case, we want to apply a bandpass filter from 0.01 to 0.08Hz

But when we look at filter response, it seems that the filter corresponds to a bandpass filter from 0.01 to 0.3Hz instead of from 0.01Hz to 0.08Hz.

whereas the version from 2016 to 2018 seems to more accurate :

Apparently, it is because of the choice of the transition band as the 'bug' only appears if we left the default transition band but disappear if we specify a correct one such as 0.001Hz

Best regards, Edouard

[sbaillet]

Hi Edouard:

Try with the 40dB relaxed option. Your bandwidth request is indeed extremely narrow and is challenging to accommodate together with sharp 60dB attenuation. I hope this will help.

Best wishes,

Sylvain.

[ftadel]

@HosseinShahabi can you please comment on the choice of the default transition band?

@Edouard2laire Do you realize that the filter you are designing here is maybe not something that you can practically apply to your data? The order of the filter is 7252: at 10Hz, this represents 725 seconds. It starts filtering correctly only after 6 minutes of recordings, and only if there are 6 extra minutes at the end of the file... I would not consider preventing this from happening to be a bug, but a decent safeguard. When you have enough information, please close the issue.

[Edouard2laire]

Hello, Thanks for your responses. I agree this is an extremely narrow bandwidth but it's very common in NIRS to filter between 0.01 to either 0.08, 0.1 or 0.15Hz and we face the same issue.

With the 2016-208 filter version, when we filter from 0.01 to 0.08, the windows say that the transient window that contains 99% of the energy is 55s for 60db attenuation or 49s for 40db. So I found this confusing with the fact that the first 6 minutes are not filtered properly. Would zero-padding the signal be able to help us in that situation?

I would not consider preventing this from happening to be a bug, but a decent safeguard. When you have enough information, please close the issue.

Maybe adding a warning, in this case, might be a good idea :)

[ftadel]

@HosseinShahabi Can you please address this?

[HosseinShahabi]

@Edouard2laire I'm checking this to see what are our default values but as you said it is a very narrow filter and usually we cannot do better than that with an FIR filter. If we are allowed to use an IIR filter, that's a different story.

Please consider that the warning is too conservative in general and in many cases it doesn't happen. I had some simulations in the past so I will post them if I can find them.

Meanwhile, would you please reference the method section of other NIRS papers so I can take a look at them and see how they did filtering?

Thanks and I will get back to you soon.

[Edouard2laire]

If we are allowed to use an IIR filter, that's a different story.

Indeed, IIR seems to be frequent (Butterworth). Following article (https://www.mdpi.com/1999-4893/11/5/67/htm) say :
Several different filters have been implemented in fNIRS studies for noise removal, with the majority using the Butterworth filter. Other filters used include Chebyshev, elliptic, and zero-phase Fast Fourier Transform (FFT) filters.

The main frequency band analyzed in fNIRS are detailed here : https://ieeexplore.ieee.org/abstract/document/8037377

Physiological fluctuations are commonly present in functional studies of hemodynamic response such as functional magnetic resonance imaging (fMRI) and functional near-infrared spectroscopy (fNIRS). These studies have revealed a distinct spectrum of these fluctuations in the frequency band of 0.01 to 2 Hz [3]–[4][5]. While the main sources of physiological fluctuations in the range of 0.15 to 2 Hz are cardiac pulsation (~ 1 Hz) and respiration (~ 0.3 Hz), low frequency oscillations in the range of 0.01 to 0.15 Hz remain controversial and these fluctuations have been attributed to myogenic (~ 0.1 Hz) and neurogenic (~ 0.02 Hz) activities [6] [7]. Evidence suggests that these physiological signals may correlate with the stimulation response in the same brain regions as the activation target which may lead to incorrect interpretation of fNIRS activation maps. Therefore, these fluctuations are commonly treated as noise and are reduced through averaging, filtering, or regression analysis to minimize their effect on fNIRS signals.

HS: Is it possible to use the ICA version of brainstorm implemented for EEG to clean NIRS data?

[ftadel]

Is it possible to use the ICA version of brainstorm implemented for EEG to clean NIRS data?

Brainstorm uses the runica.m function from EEGLAB. We are not experts in this field, you'd better ask on the EEGLAB side.

[HosseinShahabi]

@Edouard2laire Thanks for sharing those articles. We can think about using IIR filters which might have a smaller practical transient effect. We have used them in band-stop and notch filters, where the non-linear phase appears in the stopband. However, we wanna make sure that the possible distortion is minimal in that case.

Also, we already have the Morlet wavelet in Brainstorm and using a filter bank can be another approach.

For ICA, I'm a little bit skeptical that ICA can be able to decompose the signal to those specific frequency bands you are mentioning.

I need to talk with our team to see what's the best strategy to handle this issue.

[Edouard2laire]

Hi Thanks for your reply. I just found this very nice article (https://www.ncbi.nlm.nih.gov/pmc/articles/PMC6336925/) that does a review of the diverse signal filtering methods and the frequency band used in fNIRS and who indeed seems to suggest to use IIR filter. I will read it carefully, do some test and come back if I have a specific question :)

@ftadel what would be the better solution in this case: having a separate process for IIR filter or adding an option in the current process to chose between IIR and FIR filter?

[ftadel]

This is a more a question for @HosseinShahabi

[Edouard2laire]

ok. So I have made some change to the current process that exist in Nirstorm ( https://github.com/Nirstorm/nirstorm/pull/116, cc @thomas-vincent ) which is marked as deprecated.

Therefore, Here is the filter specification implented in NIRSTORM (order 3 and 4):

It has two potential problems:

• the filter is unstable after order 3. Even when using the same frequency band as the article previously linked ( [https://www.ncbi.nlm.nih.gov/pmc/articles/PMC6336925/figure/F6/] figure 6), order 4 bandpass filter isn't stable :
  • it has a large edge effect (but, I guess we just have to remove the time segment from the file after filtering?) note: The process is using filtfilt to avoid the phase issue.

Do you have any idea why the filter is unstable and how we could manage this? Is there a way to reduce the edge effect of the filter? Is it a problem if the filter already have an attenuation of 5db at the cutoff frequency?

Thanks.

[ftadel]

@HosseinShahabi @sbaillet @richardmleahy ?

[HosseinShahabi]

@Edouard2laire I'm sorry that I missed this conversation.

So, is this the filter you designed based on the Frontiers paper? https://github.com/Nirstorm/nirstorm/blob/master/bst_plugin/process_nst_iir_filter.m

I'm checking the code to see how it works.

It is very predictable that an IIR becomes unstable. The solution is to relax the transition band in this case since it is a Butterworth filter. This type of filter is more useful when you need an extremely flat response in the passband.

But if you cannot accept a relax transition band, you might better to switch to Chebyshev filter, which has a sharper transition band but accept a bit ripple in the passband.

I also checked the paper you mentioned. As you see, the majority of papers did not report the filter type (it might be Chebyshev) while Butterworth is the most common IIR one as I see.

You might also consider that some previous papers might have had similar problems as you have right now and they could not achieve the most accurate filters, which is fine in some cases.

So, in case you are looking for a sharp transition band (if I remember correctly), the Chebyshev or Elliptic might be better options.

[ftadel]

@HosseinShahabi @Edouard2laire Then what needs to be implemented in Brainstorm? What does the interface need to look like? (parameters of interest, etc) Does it it need to be a separate process, or is it an additional option you could add to the existing band-pass filter process?

[Edouard2laire]

I am sorry, I didn't had much time to work on this.

@HosseinShahabi Actually, this filter have been designed by @thomas-vincent a long time ago, I just changed recently the output and added some visualization to make it similar to the process that is in brainstorm. But I think, yes it correspond to what is in the article (a third order Butterworth filter)

From what I understood, the main characteristic we need for our filter are :

• constant gain in the bandpass as we don't want to favorise specific frequency
• constant / no phase delay

From what I understood, @ftadel, from out lab meeting, what we need to have this implanted in brainstorm is to :

• make the default transition band as a % of the bandwidth of the filter instead of a fix value so its still ok nirs

  Then you can add an expert mode that allows to change the following options :

• choice between IIR and FIR filter.
• choice to keep the mean or not even when doing a bandpass filter.

For the IIR filter design, we seems to have two options to fit the Butterworth filter:

• what have been done by @thomas-vincent here : https://github.com/Nirstorm/nirstorm/blob/master/bst_plugin/process_nst_iir_filter.m which is to specify the bandpass windows and filter order (apparently 3 in our case) and to use Butter to fit the filter to the specification, then use filtfilt to apply it to the data (it make the filter zero-pahse as I understood)

• MNE seems to specify the bandpass windows and the transition band (which seems also more like what is done in brainstorm for FIR filter) but I don't know how they fit the filter to the filter specification : https://mne.tools/mne-nirs/auto_examples/plot_20_cui.html#import-and-preprocess-data

Actually, I don't now why we are using such a narrow band here. Maybe we will need to adjust the tutorial but I need to do more test and discuss that with our lab first :)

[Edouard2laire]

Dear all :)

You were right, it was mainly our usage of the filter that was wrong. The main problem was the high-pass filter used to remove trends and slow-fluctuations. (eg a period larger to 200s or a frequency inferior to 0.05Hz) . (causing the design to be very hard to fit)

So I did the test to remove those slow frequencies outside the band-pass filter. I did some tests using brainstorm detrending which worked great and added in nirstorm an implementation to also remove discrete cosinus term to remove other slow fluctuations (eg from a period of 200s to the acquisition duration). I can send you examples if you want, but I don't want this post to be too long.

So after that, we only need the FIR to accomplish the low-pass filtering, and then using a low cutoff frequency of 0.08Hz and a transition band of 0.02 (and 40db). Then we have nice results with still a low order filter. (eg 1min transient). And we still can use the IIR lowpass filter in NIRSTORM. So I think it's not useful to spend more time on development to merge it in Brainstorm as we will be the only users? what do you think?

The only issue is that when using the FIR filter with nirs, we have to manually specify the transition band. As the default one is too large. If you can change that to a % value of the bandpass windows it's great. otherwise, I think we can just specify in the tutorial that this value have to be set.

Thanks a lot for your help finding our issue :)
You can close this issue if you want @ftadel

Edit: also, I created my account on the website to edit tutorial. (my name is EdouardDelaire) Could you allow me to create pages and upload images ? So I will update the nirs tutorial, is it better for you that I make the changes directly in the current tutorial or that I create a new one ?

Best regards, Edouard

[ftadel]

@Edouard2laire Sorry I didn't reply to this earlier. I'm closing the issue now as you requested it, but don't hesitate to reopen it you think anything from NIRSTORM could be interesting to be merged into the Brainstorm distribution.

Could you edit the tutorial page the way you wanted?