Merge pull request #128 from florisvb/notebooks-update-for-optimizer

Notebooks update for optimizer

Author: pavelkomarov

examples/2a_optimizing_parameters_with_dxdt_known.ipynb

@@ -9,3 +9,4 @@
 from pynumdiff.linear_model import lineardiff, polydiff, spectraldiff, savgoldiff
 from pynumdiff.kalman_smooth import constant_velocity, constant_acceleration, constant_jerk,\
     known_dynamics
+from pynumdiff.optimize import optimize

examples/2b_optimizing_parameters_with_dxdt_unknown.ipynb

@@ -39,7 +39,7 @@ def savgoldiff(x, dt, params=None, options=None, poly_order=None, window_size=No
         raise ValueError("`poly_order`, `window_size`, and `smoothing_win` must be given.")
 
     window_size = np.clip(window_size, poly_order + 1, len(x)-1)
-    if not window_size % 2: window_size += 1 # window_size needs to be odd
+    if window_size % 2 == 0: window_size += 1 # window_size needs to be odd
     smoothing_win = min(smoothing_win, len(x)-1)
 
     dxdt_hat = scipy.signal.savgol_filter(x, window_size, poly_order, deriv=1)/dt
@@ -312,6 +312,9 @@ def _lineardiff(x, dt, order, gamma, solver=None):
 
     if not window_size:
         return _lineardiff(x, dt, order, gamma, solver)
+    elif window_size % 2 == 0:
+        window_size += 1
+        warn("Kernel window size should be odd. Added 1 to length.")
 
     kernel = {'gaussian':utility.gaussian_kernel, 'friedrichs':utility.friedrichs_kernel}[kernel](window_size)
 

pynumdiff/__init__.py

@@ -20,11 +20,12 @@
                    'pad_to_zero_dxdt': True,
                    'high_freq_cutoff': [1e-3, 5e-2, 1e-2, 5e-2, 1e-1]},
                   {'high_freq_cutoff': (1e-5, 1-1e-5)}),
-    polydiff: ({'step_size': 1,
+    polydiff: ({'step_size': [1, 2, 5],
                 'kernel': 'friedrichs',
                 'poly_order': [2, 3, 5, 7],
                 'window_size': [11, 31, 51, 91, 131]},
-               {'poly_order': (1, 8),
+               {'step_size': (1, 100),
+                'poly_order': (1, 8),
                 'window_size': (10, 1000)}),
     savgoldiff: ({'poly_order': [2, 3, 5, 7, 9, 11, 13],
                   'window_size': [3, 10, 30, 50, 90, 130, 200, 300],
@@ -41,7 +42,8 @@
                   'order': 3,
                   'gamma': [1e-1, 1, 10, 100],
                   'window_size': [10, 30, 50, 90, 130]},
-                 {'gamma': (1e-3, 1000),
+                 {'order': (1, 5),
+                  'gamma': (1e-3, 1000),
                   'window_size': (15, 1000)}),
     first_order: ({'num_iterations': [5, 10, 30, 50]},
                   {'num_iterations': (1, 1000)}),
@@ -130,8 +132,8 @@ def optimize(func, x, dt, search_space={}, dxdt_truth=None, tvgamma=1e-2, paddin
                     structured as :code:`{param1:[values], param2:[values], param3:value, ...}`. The search space
                     is the Cartesian product. If left None, a default search space of initial values is used.
     :param np.array[float] dxdt_truth: actual time series of the derivative of x, if known
-    :param float tvgamma: regularization value used to select for parameters that yield a smooth derivative.
-                    Larger value results in a smoother derivative
+    :param float tvgamma: Only used if :code:`dxdt_truth` is given. Regularization value used to select for parameters
+                    that yield a smooth derivative. Larger value results in a smoother derivative.
     :param int padding: number of time steps to ignore at the beginning and end of the time series in the
                     optimization. Larger value causes the optimization to emphasize the accuracy of dxdt in the
                     middle of the time series

pynumdiff/linear_model/_linear_model.py

@@ -106,8 +106,8 @@ def test_spectraldiff():
     assert np.isclose(params2['high_freq_cutoff'], 0.575, atol=1e-3) or np.isclose(params2['high_freq_cutoff'], 0.735, atol=1e-3) # this one solves unstably to one or the other
 
 def test_polydiff():
-    params1, val1 = optimize(polydiff, x, dt, tvgamma=tvgamma, dxdt_truth=dxdt_truth)
-    params2, val2 = optimize(polydiff, x, dt, tvgamma=0, dxdt_truth=None)
+    params1, val1 = optimize(polydiff, x, dt, search_space={'step_size':1}, tvgamma=tvgamma, dxdt_truth=dxdt_truth)
+    params2, val2 = optimize(polydiff, x, dt, search_space={'step_size':1}, tvgamma=0, dxdt_truth=None)
     assert (params1['poly_order'], params1['window_size']) == (6, 50)
     assert (params2['poly_order'], params2['window_size']) == (4, 10)
 

pynumdiff/optimize/_optimize.py

@@ -132,6 +132,6 @@ def total_variation(x, padding=0):
     if padding == 'auto':
         padding = int(0.025*len(x))
         padding = max(padding, 1)
-    x = x[padding: len(x)-padding]
+    x = x[padding:len(x)-padding]
 
-    return np.sum(np.abs(x[1:]-x[:-1]))/len(x) # normalized version of cvxpy.tv(x2-x1).value
+    return np.linalg.norm(x[1:]-x[:-1], 1)/len(x) # normalized version of what cvxpy.tv does