Skip to content

Adds SampleSuccGauss FailureHandler for UKI#118

Merged
bors[bot] merged 1 commit intomainfrom
successful_gaussian_uki
Mar 3, 2022
Merged

Adds SampleSuccGauss FailureHandler for UKI#118
bors[bot] merged 1 commit intomainfrom
successful_gaussian_uki

Conversation

@ilopezgp
Copy link
Contributor

@ilopezgp ilopezgp commented Mar 3, 2022

Adds the Sample-successful-gaussian Failure Handler for the Unscented Kalman Inversion.

The handler is applied during the analysis step, by computing the ill-defined means and covariances over the successful particles. In order to preserve the sum of the weighting scheme, the weights at each successful sigma point are scaled up.

Consider a set of sigma points {θ_j}={θ_succ_j} U {θ_fail_j}, from which we have to estimate E[G], Cov(G, G) and Cov(G, θ) for the UI update rules. These would generally be computed through the quadrature rules of UKI (possibly modified following Huang et al 2021).

For the mean of the forward model eval we do:

  • E[G] ≈ ∑_succ [(∑W_j) / (∑_succ W_j) * W_j * G(θ_succ_j)] if the center particle did not fail (rescaled quadrature).
  • E[G] ≈ ∑_succ G(θ_succ_j) / len({θ_succ_j}) if the center particle failed (arithmetic mean).

The second option is needed because in the modified unscented transform, W_1 = 1 and all other weights are zero. We do similarly for the mean over parameters θ.

For the covariances we scale up the weights of the off-centered sigma points and compute over the successful particles,

  • Cov(G, θ) = ∑_succ [(∑_(j>1) Wc_j) / (∑_(succ, j>1) Wc_j) * Wc_j * (G(θ_succ_j) - G(m)) * (θ_succ_j - m)]
  • Cov(G, θ) = ∑_succ [(∑_(j>1) Wc_j) / (∑_(succ, j>1) Wc_j) * Wc_j * (θ_succ_j - m) * (θ_succ_j - m)]
  • Cov(G, G) = Wc_1 * (G(m) -E[ G]) * (G(m) -E[ G]) + ∑_succ_(j>1) [(∑_(j>1) Wc_j) / (∑_(succ, j>1) Wc_j) * Wc_j * (G(θ_succ_j) - E[G]) * (G(θ_succ_j) - E[G])]

Note that the term Wc_1 * (G(m) -E[ G]) * (G(m) -E[ G]) disappears when using the modified unscented transform because we take E[G]≈G(m). Also, notice that the rescaling is computed over the off-centered sigma points to remain consistent for both Cov(G, θ) and Cov(G, G).

@ilopezgp ilopezgp force-pushed the successful_gaussian_uki branch from 6af287f to 874ab6b Compare March 3, 2022 06:09
@ilopezgp ilopezgp force-pushed the successful_gaussian_uki branch from 874ab6b to 665a155 Compare March 3, 2022 17:35
Copy link
Contributor

@Zhengyu-Huang Zhengyu-Huang left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The improvement looks great! Thanks!

Copy link
Member

@odunbar odunbar left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Looks good. Two points

  • can we guarantee that the IgnoreFailures() will pass through and produce NaNs rather than cause errors? (which is what is tested here)
  • can we have some warning flags for >50% fails etc? (as in EKI) - I presume this is still important. It may be worth at some point exporting this out of the individual methods.

@ilopezgp
Copy link
Contributor Author

ilopezgp commented Mar 3, 2022

Looks good. Two points

  • can we guarantee that the IgnoreFailures() will pass through and produce NaNs rather than cause errors? (which is what is tested here)
  • can we have some warning flags for >50% fails etc? (as in EKI) - I presume this is still important. It may be worth at some point exporting this out of the individual methods.

Thanks! About the points:

  • We can guarantee that IgnoreFailures() will pass through and produce NaNs during a single call to update_ensemble!. However, when the new parameters (NaNs) are passed to the forward model, the behavior is controlled by the forward model. For toy models in Julia, it would still pass through, but for more complex code in possibly other languages or involving IO, it will likely break. This is the same behavior we had prior to adding FailureHandlers, and I think at some point we may want to shift to using an actual FailureHandler as a default.
  • The 50% flag is used by all algorithms through split_indices_by_success.

@ilopezgp
Copy link
Contributor Author

ilopezgp commented Mar 3, 2022

bors r+

@bors
Copy link
Contributor

bors bot commented Mar 3, 2022

Build succeeded:

@bors bors bot merged commit a0e4f28 into main Mar 3, 2022
@bors bors bot deleted the successful_gaussian_uki branch March 3, 2022 20:33
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

3 participants