Change of default behavior for MAXIMA_THRESHOLD_WINDOW #49
Conversation
cosimolupo
commented
Sep 26, 2022
cosimolupo
commented
- MAXIMA_THRESHOLD_WINDOW has now a default value equal to 'None', corresponding to set it equal to the entire signal duration. Also accepts float values.
- Little fix involving plot_minima function in minima.py script, now plotting the interpolated position of minima when required.
|
@cosimolupo Please look through again and update. |
| @@ -200,7 +203,7 @@ def plot_minima(asig, event, channel, maxima_threshold_window, | |||
| plot_minima(asig=time_slice(asig, args.plot_tstart, args.plot_tstop), | |||
| event=time_slice(transition_event, args.plot_tstart, args.plot_tstop), | |||
| channel=int(channel), | |||
| maxima_threshold_window = args.maxima_threshold_window, | |||
| maxima_threshold_window = min(args.maxima_threshold_window, args.plot_tstop-args.plot_tstart), | |||
There was a problem hiding this comment.
This change is redundant with the changes in lines 186-188
| if args.maxima_threshold_window is None or args.maxima_threshold_window > asig.t_stop - asig.t_start: | ||
| args.maxima_threshold_window = asig.t_stop - asig.t_start |
There was a problem hiding this comment.
args.maxima_threshold_window is a float value for a time period in seconds, whereas asig.t_stop and asig.t_start are quantity values. So, before comparing them or assigning one to the other, the quantity values need to be rescaled to seconds and converted to floats.
cosimolupo
force-pushed
the
maxima_threshold_window_PR
branch
from
961e393 to
56cf073
Compare
|
I went through the proposed changes and checked for consistency and coherence; done from my side |
| @@ -64,9 +64,11 @@ NUM_INTERPOLATION_POINTS: 0 | |||
| # minimum distance between two peaks (s) | |||
| MIN_PEAK_DISTANCE: 0.28 | |||
| # amplitude fraction to set the threshold detecting local maxima | |||
| MAXIMA_THRESHOLD_FRACTION: .5 | |||
| # range of allowed values is [0,1], default value is 0.5 | |||
There was a problem hiding this comment.
The comments shouldn't describe default values, except they have special significance.
You can just remove this line, since the range is already implied by 'fraction'
There was a problem hiding this comment.
Fixed with a new commit.
| def moving_threshold(signal, window, fraction): | ||
| # compute a dynamic threshold function through a sliding window | ||
| # on the signal array | ||
| strides = np.lib.stride_tricks.sliding_window_view(signal, window) | ||
| threshold_func = np.min(strides, axis=1) + fraction*np.ptp(strides, axis=1) | ||
| def moving_threshold(signal, maxima_threshold_window, maxima_threshold_fraction): | ||
| threshold_signal = [] | ||
|
|
||
| sampling_rate = signal.sampling_rate.rescale('Hz').magnitude | ||
| duration = float(signal.t_stop.rescale('s').magnitude - signal.t_start.rescale('s').magnitude) | ||
| if maxima_threshold_window is None or maxima_threshold_window > duration: | ||
| maxima_threshold_window = duration | ||
| window_frame = int(maxima_threshold_window*sampling_rate) | ||
|
|
||
| for channel, channel_signal in enumerate(signal.T): | ||
| if np.isnan(channel_signal).any(): | ||
| threshold_func = np.full(np.shape(signal)[0], np.nan) | ||
| else: | ||
| # compute a dynamic threshold function through a sliding window | ||
| # on the signal array | ||
| strides = np.lib.stride_tricks.sliding_window_view(channel_signal, window_frame) | ||
| threshold_func = np.min(strides, axis=1) + maxima_threshold_fraction*np.ptp(strides, axis=1) | ||
| # add elements at the beginning | ||
| threshold_func = np.append(np.ones(window_frame//2)*threshold_func[0], threshold_func) | ||
| threshold_func = np.append(threshold_func, np.ones(len(channel_signal)-len(threshold_func))*threshold_func[-1]) | ||
| threshold_signal.append(threshold_func) | ||
|
|
||
| threshold_signal = neo.AnalogSignal(np.array(threshold_signal).T, units=signal.units, sampling_rate=signal.sampling_rate) |
There was a problem hiding this comment.
Why create a separate threshold analogsignal here? It seems it requires more code and more RAM
There was a problem hiding this comment.
In the previous version of the code, the threshold was computed channel-wise in the detect_minima function, and then again in plot_minima function for each channel to be plotted. The new implementation does it only once for each channel.
Also, for a given channel to be visualized, the signal to be plotted was trimmed as [plot_tstart, plot_tstop] before passing it to plot_minima function. This implied that the moving threshold had to be computed on the signal after the trimming, with minima found and plotted consequently. This was different from what done in detect_minima function, where moving thresholds and minima were computed on the whole signal, with no trimming. This means that minima reported in plots could have been different from minima actually detected for a certain channel. Thanks to @mdenker for pointing this out during one of our discussions.
So, we decided to coherently compute moving threshold on the whole signal for each channel, only once, and then use such threshold_signal for coming minima detection and visualization. On one side, this certainly implies a larger RAM usage, but ensures a more coherent and robust computation of the moving thresholds. At the same time, larger number of code lines in moving_threshold function are almost balanced by code lines now useless and hence removed elsewhere in the code.
|
@rgutzen Please recheck the latest changes. |
