Stateful filtering fails with multi-dim array

[hergum]

Using more than one dimension in the input array makes filt fail for stateful filter.

s = rand(10,3)
df = digitalfilter(Lowpass(0.25), Butterworth(4))
f = DF2TFilter(df)
filt(f,s)

ERROR: MethodError: no method matching filt(::DF2TFilter{SecondOrderSections{Float64,Float64},Array{Float64,2}}, ::Array{Float64,2})

One-dimensional input works as expected. Two-dimensional inputs works along the first dimension, as expected, for stateless filtering filt(df,s).

Any ideas on how to make this work according to the doc?

[galenlynch]

You're right that the docs suggest that this should work. However, that doc string was written for normal filter coefficients, and not stateful filters. It's not really clear to me what the expected behavior would be when a stateful filter is applied to a matrix. For stateless filtering each column is treated as a separate signal, where the filter state is reset at the start of each column of the signal, and the filter state is discarded at the end of each column. If you were to do the same with a stateful filter, it's not clear what the filter state should be either at the start of each column, or after all the columns have been filtered (i.e. which column of the signal should the state reflect).

Depending on what you are expecting the filter state to do between columns of your signal, you could write an explicit loop accordingly. The only thing that really makes sense to me, given your decision to use a stateful filter, is that each column is a continuation of the same signal. In this case, there is no reason to use a matrix instead of a vector. If you instead want to maintain separate filter state for each column, then you should make separate filter objects for each column.

Could you explain your use case a little bit more? Are you processing very long signals in batches, or are you doing some on-line filtering of signals? Are you filtering multiple channels of data that should be treated separately? What is the reason that you want to use stateful filters in the first place?

[hergum]

I would like the stateful filter to keep track of the state of the filter for each column.

My use-case is I get long 2D signals in batches. Each row is independent, and all rows should be processed the same way. The rows are continuous in time, and must be processed as if I had the full signal available, so I need to preserve the filter state. Its for a medical device, but the use case might as well be wanting to low-pass filter the audio track of 15 radio stations. In matlab I'm used to doing things like:

>> [b,a]=butter(4,0.25);
>> [s,state]=filter(b,a,rand(15,800),[],2); % [] means no state for the first batch
>> size(state)

ans =

     4    15
>> [s,state]=filter(b,a,rand(15,800),state,2);  % use the filter state from the end of last batch to initialize the next

The documentation suggested I could do this very conveniently in Julia with a stateful filter. I'll do it in a loop instead.

[galenlynch]

I see, that makes sense. Unfortunately the design of stateful filters in this package requires that they pre-allocate the state associated with the filter. This makes them rigid, and incapable of handling either changing the type of input that the can handle (your last issue), or the number of inputs (this issue). It seems like the Matlab approach of simply returning the state along with the filtered result is simple, and also more flexible. The down side of doing something similar in Julia is that functions can't change their number of outputs as easily as in Matlab, meaning that if we were to always return the state, then everyone using `filt` would have to handle the state output. This might still be worth the trade off, both because it would simplify the code, and eliminate both of the valid issues that you've raised. That being said, I've always found `filt` treating multi-dimensional arrays as many signals, instead of only operating on vectors, to be a holdover from Matlab, and haven't been sure it's a great idea in the first place.

…

On Mon, Aug 17, 2020, 4:49 PM hergum ***@***.***> wrote: I would like the stateful filter to keep track of the state of the filter for each column. My use-case is I get long 2D signals in batches. Each row is independent, and all rows should be processed the same way. The rows are continuous in time, and must be processed as if I had the full signal available, so I need to preserve the filter state. Its for a medical device, but the use case might as well be wanting to low-pass filter the audio track of 15 radio stations. In matlab I'm used to doing things like: >> [b,a]=butter(4,0.25); >> [s,state]=filter(b,a,rand(15,800),[],2); % [] means no state for the first batch >> size(state) ans = 4 15 >> [s,state]=filter(b,a,rand(15,800),state,2); % use the filter state from the end of last batch to initialize the next The documentation suggested I could do this very conveniently in Julia with a stateful filter. I'll do it in a loop instead. — You are receiving this because you commented. Reply to this email directly, view it on GitHub <#372 (comment)>, or unsubscribe <https://github.com/notifications/unsubscribe-auth/ABUDZTYOQV6I53TKM6T3NZ3SBGJULANCNFSM4QBN2EKQ> .

[martinholters]

Adding this to the milestone for 0.8 so that we don't forget to consider any breaking changes (like adding the state as a second return value) required for fixing this and whether we can make them land in time for 0.8. If we cannot work out what the solution should be, we can always postpone to 0.9 (or any 0.8.x if we figure out a non-breaking way).

[martinholters]

From a brief glance, I think it should be possible to make filtering with DF2TFilter work column-wise if we allow its constructor to take an extra argument for the (extra) dimensionality of the input to be filtered. E.g. for the OP example, f = DF2TFilter(df,3) would be used to create a stateful filter for that input. Without the extra argument, nothing changes. And as DF2TFilter is parametric on the type of the state, it should be possible to add this in a non-breaking way.

This issue also made me look at the extra state input for the filt(::FilterCoefficients{:z}, x, si) methods, and that does not make much sense to me. The exact usage of the state seems to be an implementation detail, so passing anything other than all-zero (of the required size) is likely to yield surprises. Unless one could pass the final state of a previous call to filt, which one cannot, as it is not returned. It also also not mutated in-place, not even by filt!, so what's this extra argument good for in the present state? Also, for the use-case discussed here, the (original) input state is replicated for each column, which also seems rarely useful. Does anyone have a good use-case for passing in si with the current implementation?

EDIT: For multi-column input, it is possible to pass in multi-column state as well.

[martinholters]

Taking off the milestone, as #602 confirms this can be solved in a non-breaking way.