Final smoothing utility changes

* remove smooth_by_geoid (was only used in google_health)
* improve documentation with more usage examples
* add a few more tests

Author: dshemetov

_delphi_utils_python/delphi_utils/__init__.py

@@ -7,7 +7,7 @@
 from .archive import ArchiveDiffer, GitArchiveDiffer, S3ArchiveDiffer
 from .export import create_export_csv
 from .utils import read_params
-from .smooth import Smoother, smoothed_values_by_geo_id
+from .smooth import Smoother
 from .geomap import GeoMapper
 
 __version__ = "0.1.0"

_delphi_utils_python/delphi_utils/smooth.py

@@ -13,19 +13,16 @@
 import warnings
 
 import numpy as np
-import pandas as pd
 
 
 class Smoother:
     """
-    This is the smoothing utility class. It handles imputation and smoothing.
-    Reasonable defaults are given for all the parameters, but fine-grained
-    control is exposed.
-
-    Instantiating a smoother class specifies a smoother with a host of parameters,
-    which can then be applied to an np.ndarray with the function smooth:
-    > smoother = Smoother(smoother_name='savgol', window_length=28, gaussian_bandwidth=144)
-    > smoothed_signal = smoother.smooth(signal)
+    This is the smoothing utility class. This class holds the parameter settings for its smoother
+    methods and provides reasonable defaults. Basic usage can be found in the examples below.
+    
+    The smoother function takes numpy arrays as input, expecting the values to come from a
+    regularly-spaced time grid. NANs are ok, as long as the array does not begin with a NAN. The
+    rest of the NANs will be handled via imputation by default, though this can be turned off.
 
     Parameters
     ----------
@@ -53,8 +50,8 @@ class Smoother:
         35                       905
         42                       1303
     impute: {'savgol', 'zeros', None}
-        If 'savgol' (default), will fill nan values with a savgol fit on the largest available time 
-        window prior (up to window_length). If 'zeros', will fill nan values with zeros. 
+        If 'savgol' (default), will fill nan values with a savgol fit on the largest available time
+        window prior (up to window_length). If 'zeros', will fill nan values with zeros.
         If None, leaves the nans in place.
     minval: float or None
         The smallest value to allow in a signal. If None, there is no smallest value.
@@ -66,6 +63,33 @@ class Smoother:
     ----------
     smooth: np.ndarray
         Takes a 1D signal and returns a smoothed version. Both arrays have the same length.
+
+    Example Usage
+    -------------
+    Example 1. Apply a rolling average smoother with a window of length 10.
+    >>> smoother = Smoother(smoother_name='moving_average', window_length=10)
+    >>> smoothed_signal = smoother.smooth(signal)
+
+    Example 2. Smooth a dataframe column.
+    >>> smoother = Smoother(smoother_name="savgol")
+    >>> df[col] = pd.Series(smoother.smooth(df[col].to_numpy()))
+
+    Example 3. Apply a rolling weighted average smoother, with 95% weight on the recent 2 weeks and
+               a sharp cutoff after 4 weeks.
+    >>> smoother = Smoother(smoother_name='savgol', poly_fit_degree=0, window_length=28,
+                          gaussian_bandwidth=144)
+    >>> smoothed_signal = smoother.smooth(signal)
+
+    Example 4. Apply a local linear regression smoother (equivalent to `left_gauss_linear`), with
+               95% weight on the recent week and a sharp cutoff after 3 weeks.
+    >>> smoother = Smoother(smoother_name='savgol', poly_fit_degree=1, window_length=21,
+                          gaussian_bandwidth=36)
+    >>> smoothed_signal = smoother.smooth(signal)
+
+    Example 5. Apply the identity function (simplifies code that iterates through smoothers _and_
+               raw data).
+    >>> smoother = Smoother(smoother_name='identity')
+    >>> smoothed_signal = smoother.smooth(signal)
     """
 
     def __init__(
@@ -95,14 +119,14 @@ def __init__(
         else:
             self.coeffs = None
 
-        SMOOTHERS = {"savgol", "left_gauss_linear", "moving_average", "identity"}
+        valid_smoothers = {"savgol", "left_gauss_linear", "moving_average", "identity"}
 
-        if self.smoother_name not in SMOOTHERS:
+        if self.smoother_name not in valid_smoothers:
             raise ValueError("Invalid smoother name given.")
 
-        IMPUTE_METHODS = {"savgol", "zeros", None}
+        valid_impute_methods = {"savgol", "zeros", None}
 
-        if self.impute_method not in IMPUTE_METHODS:
+        if self.impute_method not in valid_impute_methods:
             raise ValueError("Invalid impute method given.")
 
     def smooth(self, signal):
@@ -118,6 +142,9 @@ def smooth(self, signal):
         signal_smoothed: np.ndarray
             A smoothed 1D signal.
         """
+        if len(signal) < self.window_length:
+            raise ValueError("The window_length must be smaller than the length of the signal.")
+
         signal = self.impute(signal)
 
         if self.smoother_name == "savgol":
@@ -223,9 +250,9 @@ def left_gauss_linear_smoother(self, signal):
 
     def savgol_predict(self, signal):
         """
-        Fits a polynomial through the values given by the signal and returns the value 
+        Fits a polynomial through the values given by the signal and returns the value
         of the polynomial at the right-most signal-value. More precisely, fits a polynomial
-        f(t) of degree poly_fit_degree through the points signal[0], signal[1] ..., signal[-1], 
+        f(t) of degree poly_fit_degree through the points signal[0], signal[1] ..., signal[-1],
         and returns the evaluation of the polynomial at the location of signal[-1].
 
         Parameters
@@ -250,7 +277,7 @@ def savgol_coeffs(cls, nl, nr, poly_fit_degree, gaussian_bandwidth=100):
         Solves for the Savitzky-Golay coefficients. The coefficients c_i
         give a filter so that
             y = \sum_{i=-{n_l}}^{n_r} c_i x_i
-        is the value at 0 (thus the constant term) of the polynomial fit 
+        is the value at 0 (thus the constant term) of the polynomial fit
         through the points {x_i}. The coefficients are c_i are caluclated as
             c_i =  ((A.T @ A)^(-1) @ (A.T @ e_i))_0
         where A is the design matrix of the polynomial fit and e_i is the standard
@@ -266,7 +293,7 @@ def savgol_coeffs(cls, nl, nr, poly_fit_degree, gaussian_bandwidth=100):
         poly_fit_degree: int
             The degree of the polynomial to be fit.
         gaussian_bandwidth: float or None
-            If float, performs regression with Gaussian weights whose variance is 
+            If float, performs regression with Gaussian weights whose variance is
             the gaussian_bandwidth. If None, performs unweighted regression.
 
         Returns
@@ -396,35 +423,3 @@ def savgol_impute(self, signal):
                 )
                 signal_imputed[ix] = signal_imputed[ix - self.window_length : ix] @ coeffs
         return signal_imputed
-
-
-# TODO: this needs a test, probably
-def smoothed_values_by_geo_id(
-    df: pd.DataFrame, method="savgol", **kwargs
-) -> np.ndarray:
-    """Computes a smoothed version of the variable 'val' within unique values of 'geo_id'
-
-    Currently uses a local weighted least squares, where the weights are given
-    by a Gaussian kernel.
-
-    Parameters
-    ----------
-    df: pd.DataFrame
-        A data frame with columns "geo_id", "timestamp", and "val"
-    method: {'savgol', 'left_gauss_linear', 'moving_average'}
-        A choice of window smoother to use. Check the smoother method definitions
-        for specific parameters.
-
-    Returns
-    -------
-    np.ndarray
-        A one-dimensional numpy array containing the smoothed values.
-    """
-    smoother = Smoother(method, **kwargs)
-
-    df = df.copy()
-    df["val_smooth"] = 0
-    for geo_id in df["geo_id"].unique():
-        signal = df[df["geo_id"] == geo_id]["val"].values
-        df.loc[df["geo_id"] == geo_id, "val_smooth"] = smoother.smooth(signal)
-    return df["val_smooth"].values

_delphi_utils_python/tests/test_smooth.py

@@ -1,10 +1,6 @@
 """
-This file contains a number of smoothers left gauss filter used to smooth a 1-d signal.
-Code is courtesy of Addison Hu (minor adjustments by Maria).
-
-Author: Maria Jahja
-Created: 2020-04-16
-
+Tests for the smoothing utility.
+Authors: Dmitry Shemetov, Addison Hu, Maria Jahja
 """
 import pytest
 
@@ -167,3 +163,15 @@ def test_impute(self):
         )
         with pytest.raises(ValueError):
             imputed_signal = smoother.savgol_impute(signal)
+
+        # test window_length > len(signal)
+        signal = np.arange(20)
+        smoother = Smoother(smoother_name="savgol", boundary_method="identity", window_length=30)
+        with pytest.raises(ValueError):
+            smoothed_signal = smoother.smooth(signal)
+
+        # test the boundary methods
+        signal = np.arange(20)
+        smoother = Smoother(smoother_name="savgol", poly_fit_degree=0, boundary_method="identity", window_length=10)
+        smoothed_signal = smoother.savgol_impute(signal)
+        assert np.allclose(smoothed_signal, signal)