Skip to content
Open
Show file tree
Hide file tree
Changes from all commits
Commits
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
7 changes: 6 additions & 1 deletion validphys2/src/validphys/config.py
Original file line number Diff line number Diff line change
Expand Up @@ -1886,7 +1886,12 @@ def produce_processed_metadata_group(self, processed_data_grouping, metadata_gro
def produce_group_dataset_inputs_by_metadata(self, data_input, processed_metadata_group):
"""Take the data and the processed_metadata_group key and attempt
to group the data, returns a list where each element specifies the data_input
for a single group and the group_name
for a single group and the group_name.

Ordering
--------
Not runcard order. Groups appear in order of first appearance of a member dataset.
Wihtin each group, datasets follow runcard order.
"""
res = defaultdict(list)
for dsinput in data_input:
Expand Down
87 changes: 87 additions & 0 deletions validphys2/src/validphys/tests/test_theorycovariance.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,87 @@
"""
test_theorycovariance.py

Tests for the theory covariance matrices built in
:py:mod:`validphys.theorycovariance.construction`, and in particular for the
``nnfit_theory_covmat`` production rule, which is what ``vp-setupfit`` asks for
whenever a fit runcard contains a ``theorycovmatconfig``.
"""

import numpy as np
import pytest

from reportengine.table import savetable
from validphys.api import API

# Basename the fixture looks up when reading a user covmat from disk.
USER_COVMAT_FILENAME = "test_user_covmat.csv"
POINT_PRESCRIPTION = "3 point"


@pytest.fixture(scope="module")
def user_covmat_on_disk(tmp_path_factory, thcovmat_config):
"""Write out a "user" theory covmat and return ``(directory, matrix)``.

The matrix written out is the ``3 point`` scale variation covmat for the very
same data, so that once it has been read back by ``fromfile_covmat`` the
resulting ``user_covmat`` is numerically identical to ``theory_covmat_custom``.
This exercises the real cut and index handling of ``fromfile_covmat`` and
gives the tests below an expected value with no hardcoded numbers in it.

``theory_covmat_custom`` is float32; it is promoted to float64 before being
written, because the shortest repr of a float32 only round trips to float32
and would introduce a ~1e-8 relative error.
"""
covmat = API.theory_covmat_custom(
point_prescriptions=[POINT_PRESCRIPTION], **thcovmat_config
).astype(np.float64)

# Path must be relative to cwd because ``Loader.check_vp_output_file`` rejects absolute paths.
path = tmp_path_factory.mktemp("user_covmat")
savetable(covmat, path / USER_COVMAT_FILENAME)
return path, covmat


@pytest.mark.parametrize(
("covmat_config", "expected_factor"),
[
# Scale variations *and* a user covmat.
(
{"point_prescriptions": [POINT_PRESCRIPTION], "user_covmat_path": USER_COVMAT_FILENAME},
2.0,
),
# Scale variations only.
({"point_prescriptions": [POINT_PRESCRIPTION]}, 1.0),
# User covmat only.
({"user_covmat_path": USER_COVMAT_FILENAME}, 1.0),
# User covmat only, rescaled by ``mult_factor``.
({"user_covmat_path": USER_COVMAT_FILENAME, "mult_factor": 2.5}, 2.5),
],
ids=["scalevar_and_user", "scalevar_only", "user_only", "user_only_rescaled"],
)
def test_nnfit_theory_covmat(
monkeypatch, thcovmat_config, user_covmat_on_disk, covmat_config, expected_factor
):
"""Every branch of ``produce_nnfit_theory_covmat`` must resolve and return a
covmat indexed in runcard order on both axes to match with the experimental
covariance matrix.
"""
covmat_dir, scalevar_covmat = user_covmat_on_disk
monkeypatch.chdir(covmat_dir)
covmat = API.nnfit_theory_covmat(**covmat_config, **thcovmat_config)

runcard_index = API.data_index(**thcovmat_config).droplevel(0)
process_index = API.procs_index(**thcovmat_config).droplevel(0)

# This test only has an effect if process grouping actually reorders the runcard.
assert not runcard_index.equals(process_index), (
"Grouping DATA_THCOVMAT by process no longer reorders it, so the "
"assertions below would hold even if the covmat were written out in "
"process order instead of runcard order."
)

assert covmat.index.droplevel(0).equals(runcard_index)
assert covmat.columns.equals(covmat.index)

expected = scalevar_covmat.reindex(index=covmat.index, columns=covmat.columns)
np.testing.assert_allclose(covmat.to_numpy(), expected_factor * expected.to_numpy(), rtol=1e-8)
17 changes: 0 additions & 17 deletions validphys2/src/validphys/theorycovariance/construction.py
Original file line number Diff line number Diff line change
Expand Up @@ -17,7 +17,6 @@
from validphys.theorycovariance.higher_twist_functions import compute_deltas_pc
from validphys.theorycovariance.theorycovarianceutils import (
check_correct_theory_combination,
check_fit_dataset_order_matches_grouped,
process_lookup,
)

Expand Down Expand Up @@ -543,7 +542,6 @@ def user_covmat(


@table
@check_fit_dataset_order_matches_grouped
def total_theory_covmat(theory_covmat_custom, user_covmat):
"""
Sum of scale variation and user covmat, where both are used.
Expand Down Expand Up @@ -595,21 +593,6 @@ def user_covmat_fitting(user_covmat, data_input_matched_procs_index):
return _reindex_covmat_to_fitting_order(user_covmat, data_input_matched_procs_index)


def procs_index_matched(groups_index, procs_index):
"""procs_index but matched to the dataset order given
by groups_index."""
# Making list with exps ordered like in groups_index
groups_ds_order = groups_index.get_level_values(level=1).unique().tolist()
# Tuples to make multiindex, ordered like in groups_index
tups = []
for ds in groups_ds_order:
for orig in procs_index:
if orig[1] == ds:
tups.append(orig)

return pd.MultiIndex.from_tuples(tups, names=("process", "dataset", "id"))


@table
def theory_corrmat_custom(theory_covmat_custom):
"""Calculates the theory correlation matrix for scale variations
Expand Down
24 changes: 0 additions & 24 deletions validphys2/src/validphys/theorycovariance/theorycovarianceutils.py
Original file line number Diff line number Diff line change
Expand Up @@ -85,30 +85,6 @@ def check_correct_theory_combination_internal(
check_correct_theory_combination = make_argcheck(check_correct_theory_combination_internal)


@make_argcheck
def check_fit_dataset_order_matches_grouped(
group_dataset_inputs_by_metadata, data_input, processed_metadata_group
):
"""
Check for use with theory covmat generation.

Makes sure that the order of datasets listed in the fit runcard is the same
as that specified by the metadata grouping. Otherwise there can be a
misalignment between the experiment covmat and theory covmat.
"""
data_input_iter = iter(data_input)
for group in group_dataset_inputs_by_metadata:
for dsinput in group["data_input"]:
grouped_ds = dsinput.name
input_ds = next(data_input_iter).name
check(
grouped_ds == input_ds,
"Dataset ordering is changed by grouping, this will cause "
"errors when running fits with theory covmat. Datasets should "
f"be ordered by {processed_metadata_group} in the runcard.",
)

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

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

I see this check is from the time the fitting code was c++ so I'm ready to believe it is obsolete.
However, I also know that two different PhD students three year apart spent many hours on a misordered sum of covmats.
Now it would be about 3 years later from last time... are you absolutely sure this is safe to remove?

Btw, a way to check whether the hypothesis is correct from what you say would be to simply add groups_index to data_input_matched_procs_index as a dummy argument. That should be enough to trigger the production.

Copy link
Copy Markdown
Contributor Author

Choose a reason for hiding this comment

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

Btw, a way to check whether the hypothesis is correct from what you say would be to simply add groups_index to data_input_matched_procs_index as a dummy argument. That should be enough to trigger the production.

Yes, I tried to add groups_index as a dummy argument and it worked. However, I'm not in favour of this solution for two reasons. First, it's very hard to understand the reason why the dummy variable has been included. However, a comment may help in this case. Second, as you said we'd have an important part of the code (checking ordering of covmats) that rely on obselete code. Moreover, if we were to keep this obselete code, we'd have two functions checking the same thing in two not-so-different branches of the code.

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

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

No, I don't think it should be added as an argument of course.
But it is a way to check whether it was a possible issue or the actual issue.

Copy link
Copy Markdown
Contributor Author

Choose a reason for hiding this comment

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

Good. I reach these conclusions:

  • check_fit_dataset_order_matches_grouped checks if group_dataset_inputs_by_metadata follows or not runcard order.
  • However, group_dataset_inputs_by_metadata is never used to construct an index in the theory covmat pipeline.
  • The method itself is called in
    def produce_group_dataset_inputs_by_process(self, data_input):
    return self.produce_group_dataset_inputs_by_metadata(data_input, "nnpdf31_process")

    but it's a different thing.
  • Anyway, this checked was a guard and didn't fix the order of the code: it just refused to run. Now the code fixes the ordering by reindexing. On the writing side, theory_covmat_custom_fitting / total_theory_covmat_fitting / user_covmat_fitting push the covmat through _reindex_covmat_to_fitting_order with data_input_matched_procs_index. This latter matches with the order of the runcard, which is the same as the exp. covmat. On the reading side, produce_loaded_theory_covmat always reindexes.
  • The critical part is
    @table
    def total_theory_covmat(theory_covmat_custom, user_covmat):
    """
    Sum of scale variation and user covmat, where both are used.
    """
    return theory_covmat_custom + user_covmat

    where both matrices are added. However, fromfile_covmat uses the same index as in theory_covmat_custom_per_prescription, namely procs_index. If their index mismatches, then the sum would complain rasing an error.

Given these points, I think it's safe to remove this function.


def process_lookup(name):
"""
Returns the `nnpdf31_process` of the corresponding dataset.
Expand Down
Loading