Skip to content

gh-293: ensure delete-2 JK remains posdef #292

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.

Sign up for GitHub

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

Open
JaimeRZP wants to merge 2 commits into
base: main
Choose a base branch
from
Open

gh-293: ensure delete-2 JK remains posdef #292

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
aeba96c
bug
JaimeRZP Jul 4, 2025
6c2a417
ruff
JaimeRZP Jul 4, 2025
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
29 changes: 12 additions & 17 deletions heracles/dices/jackknife.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -327,15 +327,6 @@ def delete2_correction(cls0, cls1, cls2):
Q_ii.append(qii)
# Compute the correction from the ensemble
Q = _jackknife_covariance(Q_ii, nd=2)
# Diagonalise the correction
for key in Q:
q = Q[key]
*_, length = q.shape
q_diag = np.diagonal(q, axis1=-2, axis2=-1)
q_diag_exp = np.zeros_like(q)
diag_indices = np.arange(length) # Indices for the diagonal
q_diag_exp[..., diag_indices, diag_indices] = q_diag
Q[key] = Result(q_diag_exp, axis=q.axis, ell=q.ell)
Comment on lines -330 to -338
Copy link
Contributor

@ntessore ntessore Jul 7, 2025

Choose a reason for hiding this comment

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

Why is this no longer necessary here? Doesn't the change also change the result of this function?

Copy link
Contributor Author

@JaimeRZP JaimeRZP Jul 8, 2025

Choose a reason for hiding this comment

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

Indeed it does but I think this way it is better. Rather than diagonalizing the correction and then simply adding it to the cov_jk, we can save the fully correction (in case somebody wishes to have a look) and then substract its diagonal to the diagonal of the cov_jk.

return Q


Expand All @@ -359,12 +350,16 @@ def _debias_covariance(cov_jk, Q):
Internal method to debias the Jackknife covariance.
Useful when the delete2 correction is already computed.
"""
debiased_cov = {}
for key in list(cov_jk.keys()):
c = cov_jk[key].array - Q[key].array
debiased_cov[key] = Result(
c,
ell=cov_jk[key].ell,
axis=cov_jk[key].axis,
)
debiased_cov = deepcopy(cov_jk)
Copy link
Contributor

@ntessore ntessore Jul 7, 2025

Choose a reason for hiding this comment

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

I don't think the deepcopy is useful here, if you set all entries to a new Result down below anyway.

Copy link
Contributor Author

@JaimeRZP JaimeRZP Jul 8, 2025

Choose a reason for hiding this comment

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

I think otherwise the entries of cov_jk get accidentally overwritten in the process.

for key in Q:
q = Q[key].array
c = debiased_cov[key].array
*_, length = q.shape
q_diag = np.diagonal(q, axis1=-2, axis2=-1)
c_diag = np.diagonal(c, axis1=-2, axis2=-1)
# Indices for the diagonal
diag_indices = np.arange(length)
# abs is only needed when too few Jackknife regions are used
c[..., diag_indices, diag_indices] = abs(c_diag - q_diag)
Comment on lines +362 to +363
Copy link
Contributor

@ntessore ntessore Jul 7, 2025

Choose a reason for hiding this comment

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

Isn't this materially changing the output of the debiasing step? In the sense that the absolute value is something we are making up on the spot to fix a "problem" that is really intrinsic to the process?

Copy link
Contributor Author

@JaimeRZP JaimeRZP Jul 8, 2025

Choose a reason for hiding this comment

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

for large enough jk regions, q_diag should always be smaller than c_diag. The goal here is to cover for when only a few regions are used. We could something more kosher such as enfocing pos-def for each entry in the matrix stack but this would require more involved coding.

debiased_cov[key] = Result(c, axis=cov_jk[key].axis, ell=cov_jk[key].ell)
return debiased_cov
14 changes: 12 additions & 2 deletions tests/test_dices.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -2,6 +2,7 @@
import numpy as np
import heracles
import pytest
from copy import deepcopy
import heracles.dices as dices
from heracles.healpy import HealpixMapper
from heracles.fields import Positions, Shears, Visibility, Weights
Expand Down Expand Up @@ -313,9 +314,18 @@ def test_debiasing(data_path):
cqs1,
cqs2,
)
_debiased_cov = {}
_debiased_cov = deepcopy(cov_jk)
for key in list(cov_jk.keys()):
_debiased_cov[key] = cov_jk[key].array - Q[key]
q = Q[key].array
c = _debiased_cov[key].array
*_, length = q.shape
q_diag = np.diagonal(q, axis1=-2, axis2=-1)
c_diag = np.diagonal(c, axis1=-2, axis2=-1)
# Indices for the diagonal
diag_indices = np.arange(length)
# abs is only needed when too few Jackknife regions are used
c[..., diag_indices, diag_indices] = abs(c_diag - q_diag)
_debiased_cov[key] = c
Comment on lines -318 to +328
Copy link
Contributor

@ntessore ntessore Jul 7, 2025

Choose a reason for hiding this comment

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

Isn't this test the exact same code used in the function? If so, can we have a test where you construct a known output for a known input by hand? Otherwise, we might just end up with the same bug in 2 places...

Copy link
Contributor Author

@JaimeRZP JaimeRZP Jul 8, 2025

Choose a reason for hiding this comment

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

So was before to be honest


# Check diagonal
for key in list(debiased_cov.keys()):
Expand Down