Skip to content

Commit 599d5a4

Browse files
committed
deduplicated unused rmse function in evaluate
1 parent a1d7aa7 commit 599d5a4

3 files changed

Lines changed: 60 additions & 133 deletions

File tree

pynumdiff/optimize/__init__.py

Lines changed: 0 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -27,7 +27,3 @@
2727
from pynumdiff.optimize import total_variation_regularization
2828
from pynumdiff.optimize import linear_model
2929
from pynumdiff.optimize import kalman_smooth
30-
31-
32-
33-

pynumdiff/utils/evaluate.py

Lines changed: 54 additions & 121 deletions
Original file line numberDiff line numberDiff line change
@@ -1,9 +1,9 @@
11
"""
22
Metrics and evaluations?
33
"""
4-
import numpy as _np
5-
import matplotlib.pyplot as _plt
6-
import scipy.stats as _scipy_stats
4+
import numpy as np
5+
import matplotlib.pyplot as plt
6+
from scipy import stats
77

88
# local imports
99
from pynumdiff.utils import utility as _utility
@@ -13,49 +13,26 @@
1313
# pylint: disable-msg=too-many-locals, too-many-arguments
1414
def plot(x, dt, x_hat, dxdt_hat, x_truth, dxdt_truth, xlim=None, ax_x=None, ax_dxdt=None,
1515
show_error=True, markersize=5):
16-
"""
17-
Make comparison plots of 'x (blue) vs x_truth (black) vs x_hat (red)' and
16+
"""Make comparison plots of 'x (blue) vs x_truth (black) vs x_hat (red)' and
1817
'dxdt_truth (black) vs dxdt_hat (red)'
1918
20-
:param x: array of noisy time series
21-
:type x: np.array (float)
22-
23-
:param dt: a float number representing the step size
24-
:type dt: float
25-
26-
:param x_hat: array of smoothed estimation of x
27-
:type x_hat: np.array (float)
28-
29-
:param dxdt_hat: array of estimated derivative
30-
:type dxdt_hat: np.array (float)
31-
32-
:param x_truth: array of noise-free time series
33-
:type x_truth: np.array (float)
34-
35-
:param dxdt_truth: array of true derivative
36-
:type dxdt_truth: np.array (float)
37-
38-
:param xlim: a list specifying range of x
39-
:type xlim: list (2 integers), optional
40-
41-
:param ax_x: axis of the first plot
42-
:type ax_x: :class:`matplotlib.axes`, optional
43-
44-
:param ax_dxdt: axis of the second plot
45-
:type ax_dxdt: :class:`matplotlib.axes`, optional
46-
47-
:param show_error: whether to show the rmse
48-
:type show_error: boolean, optional
49-
50-
:param markersize: marker size of noisy observations
51-
:type markersize: int, optional
19+
:param np.array[float] x: array of noisy data
20+
:param float dt: a float number representing the step size
21+
:param np.array[float] x_hat: array of smoothed estimation of x
22+
:param np.array[float] dxdt_hat: array of estimated derivative
23+
:param np.array[float] x_truth: array of noise-free time series
24+
:param np.array[float] dxdt_truth: array of true derivative
25+
:param list[int] xlim: a list specifying range of x
26+
:param matplotlib.axes ax_x: axis of the first plot
27+
:param matplotlib.axes ax_dxdt: axis of the second plot
28+
:param bool show_error: whether to show the rmse
29+
:param int markersize: marker size of noisy observations
5230
5331
:return: Display two plots
54-
:rtype: None
5532
"""
56-
t = _np.arange(0, dt*len(x), dt)
33+
t = np.arange(0, dt*len(x), dt)
5734
if ax_x is None and ax_dxdt is None:
58-
fig = _plt.figure(figsize=(20, 6))
35+
fig = plt.figure(figsize=(20, 6))
5936
ax_x = fig.add_subplot(121)
6037
ax_dxdt = fig.add_subplot(122)
6138

@@ -89,125 +66,81 @@ def plot(x, dt, x_hat, dxdt_hat, x_truth, dxdt_truth, xlim=None, ax_x=None, ax_d
8966
print('RMS error in velocity: ', rms_dxdt)
9067

9168

92-
def __rms_error__(a, e):
93-
"""
94-
Calculate rms error
69+
def _rms_error(a, e):
70+
"""Calculate rms error
9571
9672
:param a: the first array
9773
:param e: the second array
98-
:return: a float number representing the rms error
74+
75+
:return: (float) -- the root mean squared error
9976
"""
100-
if _np.max(_np.abs(a-e)) > 1e16:
77+
if np.max(np.abs(a-e)) > 1e16:
10178
return 1e16
102-
s_error = _np.ravel((a - e))**2
103-
ms_error = _np.mean(s_error)
104-
rms_error = _np.sqrt(ms_error)
79+
s_error = np.ravel((a - e))**2
80+
ms_error = np.mean(s_error)
81+
rms_error = np.sqrt(ms_error)
10582
return rms_error
10683

10784

10885
def metrics(x, dt, x_hat, dxdt_hat, x_truth=None, dxdt_truth=None, padding=None):
109-
"""
110-
Evaluate x_hat based on various metrics, depending on whether dxdt_truth and x_truth are known or not.
111-
112-
:param x: time series that was differentiated
113-
:type x: np.array
114-
115-
:param dt: time step in seconds
116-
:type dt: float
117-
118-
:param x_hat: estimated (smoothed) x
119-
:type x_hat: np.array
120-
121-
:param dxdt_hat: estimated xdot
122-
:type dxdt_hat: np.array
123-
124-
:param x_truth: true value of x, if known, optional
125-
:type x_truth: np.array like x or None
126-
127-
:param dxdt_truth: true value of dxdt, if known, optional
128-
:type dxdt_truth: np.array like x or None
129-
130-
:param padding: number of snapshots on either side of the array to ignore when calculating the metric. If autor or None, defaults to 2.5% of the size of x
131-
:type padding: int, None, or auto
132-
133-
:return: a tuple containing the following:
134-
- rms_rec_x: RMS error between the integral of dxdt_hat and x
135-
- rms_x: RMS error between x_hat and x_truth, returns None if x_truth is None
136-
- rms_dxdt: RMS error between dxdt_hat and dxdt_truth, returns None if dxdt_hat is None
137-
:rtype: tuple -> (float, float, float)
138-
86+
"""Evaluate x_hat based on various metrics, depending on whether dxdt_truth and x_truth are known or not.
87+
88+
:param np.array[float] x: data that was differentiated
89+
:param float dt: step size
90+
:param np.array[float] x_hat: estimated (smoothed) x
91+
:param np.array[float] dxdt_hat: estimated xdot
92+
:param np.array[float] x_truth: true value of x, if known
93+
:param np.array[float] dxdt_truth: true value of dxdt, if known, optional
94+
:param int padding: number of snapshots on either side of the array to ignore when calculating
95+
the metric. If :code:`'auto'` or :code:`None`, defaults to 2.5% of the size of x
96+
97+
:return: tuple[float, float, float] containing\n
98+
- **rms_rec_x** -- RMS error between the integral of dxdt_hat and x
99+
- **rms_x** -- RMS error between x_hat and x_truth, returns None if x_truth is None
100+
- **rms_dxdt** -- RMS error between dxdt_hat and dxdt_truth, returns None if dxdt_hat is None
139101
"""
140102
if padding is None or padding == 'auto':
141103
padding = int(0.025*len(x))
142104
padding = max(padding, 1)
143-
if _np.isnan(x_hat).any():
144-
return _np.nan, _np.nan, _np.nan
105+
if np.isnan(x_hat).any():
106+
return np.nan, np.nan, np.nan
145107

146108
# RMS dxdt
147109
if dxdt_truth is not None:
148-
rms_dxdt = __rms_error__(dxdt_hat[padding:-padding], dxdt_truth[padding:-padding])
110+
rms_dxdt = _rms_error(dxdt_hat[padding:-padding], dxdt_truth[padding:-padding])
149111
else:
150112
rms_dxdt = None
151113

152114
# RMS x
153115
if x_truth is not None:
154-
rms_x = __rms_error__(x_hat[padding:-padding], x_truth[padding:-padding])
116+
rms_x = _rms_error(x_hat[padding:-padding], x_truth[padding:-padding])
155117
else:
156118
rms_x = None
157119

158120
# RMS reconstructed x
159121
rec_x = _utility.integrate_dxdt_hat(dxdt_hat, dt)
160122
x0 = _utility.estimate_initial_condition(x, rec_x)
161123
rec_x = rec_x + x0
162-
rms_rec_x = __rms_error__(rec_x[padding:-padding], x[padding:-padding])
124+
rms_rec_x = _rms_error(rec_x[padding:-padding], x[padding:-padding])
163125

164126
return rms_rec_x, rms_x, rms_dxdt
165127

166128

167129
def error_correlation(dxdt_hat, dxdt_truth, padding=None):
168-
"""
169-
Calculate the error correlation (pearsons correlation coefficient) between the estimated dxdt and true dxdt
130+
"""Calculate the error correlation (pearsons correlation coefficient) between the estimated
131+
dxdt and true dxdt
170132
171-
:param dxdt_hat: estimated xdot
172-
:type dxdt_hat: np.array
173-
174-
:param dxdt_truth: true value of dxdt, if known, optional
175-
:type dxdt_truth: np.array like x or None
176-
177-
:param padding: number of snapshots on either side of the array to ignore when calculating the metric. If autor or None, defaults to 2.5% of the size of x
178-
:type padding: int, None, or auto
179-
180-
:return: r-squared correlation coefficient
181-
:rtype: float
133+
:param np.array[float] dxdt_hat: estimated xdot
134+
:param np.array[float] dxdt_truth: true value of dxdt, if known, optional
135+
:param int padding: number of snapshots on either side of the array to ignore when calculating
136+
the metric. If auto or None, defaults to 2.5% of the size of x
182137
138+
:return: (float) -- r-squared correlation coefficient
183139
"""
184140
if padding is None or padding == 'auto':
185141
padding = int(0.025*len(dxdt_hat))
186142
padding = max(padding, 1)
187143
errors = (dxdt_hat[padding:-padding] - dxdt_truth[padding:-padding])
188-
r = _scipy_stats.linregress(dxdt_truth[padding:-padding] -
189-
_np.mean(dxdt_truth[padding:-padding]), errors)
144+
r = stats.linregress(dxdt_truth[padding:-padding] -
145+
np.mean(dxdt_truth[padding:-padding]), errors)
190146
return r.rvalue**2
191-
192-
193-
def rmse(dxdt_hat, dxdt_truth, padding=None):
194-
"""
195-
Calculate the Root Mean Squared Error between the estimated dxdt and true dxdt
196-
197-
:param dxdt_hat: estimated xdot
198-
:type dxdt_hat: np.array
199-
200-
:param dxdt_truth: true value of dxdt, if known, optional
201-
:type dxdt_truth: np.array like x or None
202-
203-
:param padding: number of snapshots on either side of the array to ignore when calculating the metric. If autor or None, defaults to 2.5% of the size of x
204-
:type padding: int, None, or auto
205-
206-
:return: Root Mean Squared Error
207-
:rtype: float
208-
"""
209-
if padding is None or padding == 'auto':
210-
padding = int(0.025*len(dxdt_hat))
211-
padding = max(padding, 1)
212-
RMSE = _np.sqrt(_np.mean((dxdt_hat[padding:-padding] - dxdt_truth[padding:-padding])**2))
213-
return RMSE

pynumdiff/utils/utility.py

Lines changed: 6 additions & 8 deletions
Original file line numberDiff line numberDiff line change
@@ -32,14 +32,10 @@ def hankel_matrix(x, num_delays, pad=False): # fixed delay step of 1
3232
def matrix_inv(X, max_sigma=1e-16):
3333
"""Stable (pseudo) matrix inversion using singular value decomposition. Unused throughout the repo.
3434
35-
:param X: matrix to invert
36-
:type X: np.matrix or np.array
35+
:param np.array X: matrix to invert
36+
:param float max_sigma: smallest singular values to take into account. matrix will be truncated prior to inversion based on this value.
3737
38-
:param max_sigma: smallest singular values to take into account. matrix will be truncated prior to inversion based on this value.
39-
:type max_sigma: float
40-
41-
:return: matrix pseudo inverse
42-
:rtype: np.array or np.matrix
38+
:return: (np.array) -- pseudo inverse
4339
"""
4440
U, Sigma, V = np.linalg.svd(X, full_matrices=False)
4541
Sigma_inv = Sigma**-1
@@ -51,7 +47,8 @@ def total_variation(x):
5147
"""Calculate the total variation of an array. Used by optimizer.
5248
5349
:param np.array[float] x: data
54-
:return: float total variation
50+
51+
:return: (float) -- total variation
5552
"""
5653
if np.isnan(x).any():
5754
return np.nan
@@ -162,6 +159,7 @@ def convolutional_smoother(x, kernel, iterations=1):
162159
:param np.array[float] x: 1D data
163160
:param np.array[float] kernel: kernel to use in convolution
164161
:param int iterations: number of iterations, >=1
162+
165163
:return: **x_hat** (np.array[float]) -- smoothed x
166164
"""
167165
x_hat = np.hstack((x[::-1], x, x[::-1])) # pad

0 commit comments

Comments
 (0)