Patch Smoothing

This recipe compares several smoothing strategies.

A few patch methods useful for smoothing are:

Note

Patch.rolling is quite flexible and can be used for many different processing tasks. However, as mentioned in the rolling section of the tutorial, Patch.rolling includes NaN entries in the output due to edge effects. This can require some additional thought to properly deal with.

The other methods mentioned above deal with edge effects differently (governed by the mode parameter) and by default don’t include NaN entries.

Example Data

We will use example_event_2 to demonstrate how to use each smoothing method.

import numpy as np

import dascore as dc

patch = dc.get_example_patch("example_event_2")
ax = patch.viz.waterfall()
ax.set_title("un-smoothed patch");

Rolling

Patch.rolling can be used for applying aggregation functions to rolling windows along one dimension of the Patch instance. For instance, to apply a rolling mean to the time axis:

smoothed_patch = (
    patch.rolling(time=0.005).mean()
)
ax = smoothed_patch.viz.waterfall()
ax.set_title("rolling time mean");

or a median filter to the distance axis:

smoothed_patch = (
    patch.rolling(distance=20, samples=True).median()
)
ax = smoothed_patch.viz.waterfall()
ax.set_title("rolling distance median");

Or both:

smoothed_patch = (
    patch.rolling(distance=20, samples=True).median()
    .rolling(time=0.005).mean()
)
ax = smoothed_patch.viz.waterfall()
ax.set_title("rolling time mean distance median");

Again though, the output array has NaN entries along some of the edges:

nan_data = np.isnan(smoothed_patch.data)
print(f"Number of NaN = {nan_data.sum()}")

# Create and plot a patch to highlight NaN values.
nan_patch = (
    smoothed_patch.update(data=nan_data.astype(np.int32))
    .update_attrs(data_type='', data_units=None)
)
ax = nan_patch.viz.waterfall()
ax.set_title("NaNs in Patch");
Number of NaN = 47537

Which can be dropped with Patch.dropna.

[’Patch.fill_nan`]

Savgol filter

Patch.savgol_filter uses SciPy’s savgol_filter, which is an implementation of the Savitzky-Golay Filter to apply smoothing along a single dimension.

smoothed_patch = (
    patch.savgol_filter(time=0.01, polyorder=3)
)
ax = smoothed_patch.viz.waterfall()
ax.set_title("savgol time");

Multidimensional smoothing is achieved by applying the 1D filter sequentially along each desired dimension.

smoothed_patch = (
    patch.savgol_filter(
        time=25, 
        distance=25, 
        samples=True,
        polyorder=2,
    )
)
ax = smoothed_patch.viz.waterfall()
ax.set_title("savgol time and distance");

Gaussian filter

Patch.gaussian_filter uses SciPy’s gaussian_filter to apply a gaussian smoothing operator to the patch. Note that the keyword arguments here specify standard deviation, and the truncate keyword determines how many standard deviations are included in the smoothing kernel.

smoothed_patch = patch.gaussian_filter(
    time=5, 
    distance=5, 
    samples=True,
)
ax = smoothed_patch.viz.waterfall()
ax.set_title("gaussian time and distance");

To visualize the smoothing kernel we can apply the operator to a patch whose data are all 0s except for a single 1 in the center.

data = np.zeros_like(smoothed_patch.data)
data[data.shape[0]//2, data.shape[1]//2] = 1.0

delta_patch = patch.update(data=data)

smoothed_patch = delta_patch.gaussian_filter(
    time=.005, 
    distance=15,
) 

ax = smoothed_patch.viz.waterfall()
ax.set_title("gaussian smoothing kernel");