Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Presumably shifted uncertainty result for FIRuncFilter for non-stationary uncertainty #183

Closed
mgrub opened this issue Oct 28, 2020 · 6 comments
Assignees
Labels

Comments

@mgrub
Copy link
Contributor

mgrub commented Oct 28, 2020

In cases of non-stationary uncertainty (king="diag"), the output of FIRuncFilter returns an uncertainty-result, that seems to be shifted.

Looking into the code, this is very likely caused by the equations to calculate UncCov: L.182 + L.186 + L.190 + L.192. Here, theta and Ulow are multiplied. However, for values in theta it holds theta[i+1] refers to an earlier point in time than theta[i] - while in Ulow it is reversed. Therefore, during the multiplication values are combined, that do not correspond in their time-order. To fix this, theta probably needs to be reversed, to achieve correct (time-ascending) order.

It seems like this misbehaviour is not covered by tests or examples so far. Therefore a suitable test/example should be introduced (Maybe against some MC-method). A comparison with the upcoming IIRuncFilter could also be an option.

Environment

  • OS: Windows 10
  • PyDynamic Version 1.6.0 (direct from master)

MWE showing comparison with MC-result

import matplotlib.pyplot as plt
import numpy as np

from PyDynamic.uncertainty.propagate_MonteCarlo import MC
from PyDynamic.uncertainty.propagate_filter import FIRuncFilter

n_signal = 100
n_filter = 50

x = np.append(np.arange(n_signal//2), np.zeros(n_signal//2))
ux = 0.1 * np.abs(x)

theta = np.array([0.5,0.5] + [0] * (n_filter - 2))
utheta = np.zeros((n_filter, n_filter))


# run filter
y, uy = FIRuncFilter(
    x, ux, theta, utheta, kind="diag"
)  # apply uncertain FIR filter (GUM formula)
yMC, uyMC = MC(
    x, ux, theta, [1.0], utheta, runs=1000
)  # apply uncertain FIR filter (Monte Carlo)


# compare FIR and MC results
fig, ax = plt.subplots(nrows=2)

ax[0].plot(x, label="x")
ax[0].plot(theta, label="theta")
ax[0].plot(y, label="y")
ax[0].plot(yMC, label="y (MC)")
ax[0].legend()

ax[1].plot(ux, label="ux")
ax[1].plot(np.diag(utheta), label="utheta")
ax[1].plot(uy, label="uy")
ax[1].plot(np.sqrt(np.diag(uyMC)), label="uy (MC)")
ax[1].legend()

plt.show()

Output of MWE
Current situation:
issue_on

Expected result:
issue_off

Note
The effect becomes especially obvious for long filter lengths and if the signal gets appended by many zeros.

@mgrub mgrub added the bug label Oct 28, 2020
@mgrub mgrub self-assigned this Oct 28, 2020
@mgrub mgrub changed the title Presumably reversed uncertainty result for FIRuncFilter for non-stationary uncertainty Presumably shifted uncertainty result for FIRuncFilter for non-stationary uncertainty Oct 28, 2020
@eichstaedtPTB
Copy link
Collaborator

Thanks for pointing this out. Could you extend your example above with an application of IIRuncFilter to see how that behaves in that situation?

@BjoernLudwigPTB
Copy link
Member

One first thing while starting to code on that should be implementing a test, that reveals that erroneous behaviour.

@mgrub
Copy link
Contributor Author

mgrub commented Oct 28, 2020

Thanks for pointing this out. Could you extend your example above with an application of IIRuncFilter to see how that behaves in that situation?

I will do that. However, I am not fully certain, what's the best way to achieve this. The revised IIRuncFilter currently only resides in the prepare_release_2.0.0 branch, which I don't want to pull in its entirety into a fix. Also,I can't get IIRuncFilter in the master-branch to work. So I could do the comparison in a q'n'd-fashion or we could decide to release the revised IIRuncFilter prior to PyDynamic 2.0.0.

One first thing while starting to code on that should be implementing a test, that reveals that erroneous behaviour.

Yes! However, as long as the IIRuncFilter is not available in the master branch, I can only think of an approximative check with a Monte-Carlo run. (Similar to the MWE proposed above.)

@mgrub
Copy link
Contributor Author

mgrub commented Oct 28, 2020

Here are some screenshots with the revised IIRuncFilter added to the comparison.

Issue active: IIR and MC show good agreement, FIR differs
issue_on

Issue resolved: IIR and FIR match exactly, MC shows good agreement to both
issue_resolved

Note: FIRuncFilter is based on a method, that incorporates an additional term over the pure GUM-propagation, while IIRuncFilter sticks to the GUM-scheme. Therefore IIRuncFilter and FIRuncFilter can only show the exact same output if Utheta=np.zeros((n_filter,n_filter)). As soon as Utheta becomes a full covariance matrix, the observed shift of this bug-report becomes less noticable, as the contribution of other terms becomes more prominent.

@BjoernLudwigPTB
Copy link
Member

Thanks for pointing this out. Could you extend your example above with an application of IIRuncFilter to see how that behaves in that situation?

I will do that. However, I am not fully certain, what's the best way to achieve this. The revised IIRuncFilter currently only resides in the prepare_release_2.0.0 branch, which I don't want to pull in its entirety into a fix. Also,I can't get IIRuncFilter in the master-branch to work. So I could do the comparison in a q'n'd-fashion or we could decide to release the revised IIRuncFilter prior to PyDynamic 2.0.0.

One first thing while starting to code on that should be implementing a test, that reveals that erroneous behaviour.

Yes! However, as long as the IIRuncFilter is not available in the master branch, I can only think of an approximative check with a Monte-Carlo run. (Similar to the MWE proposed above.)

Then we should simply push towards releasing as soon as possible. Or do we need a backport of this bug fix to 1.x anyway?

@mgrub
Copy link
Contributor Author

mgrub commented Oct 28, 2020

Then we should simply push towards releasing as soon as possible. Or do we need a backport of this bug fix to 1.x anyway?

I think, now that I have implemented a fix in PR #184 , we should release that in 1.x. But we should rewrite some tests after the initial 2.0.0 release.

@mgrub mgrub closed this as completed in dd04eea Oct 29, 2020
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
Projects
None yet
Development

No branches or pull requests

3 participants