diff --git a/validphys2/src/validphys/config.py b/validphys2/src/validphys/config.py index b84900de08..113e46af4b 100644 --- a/validphys2/src/validphys/config.py +++ b/validphys2/src/validphys/config.py @@ -1127,14 +1127,13 @@ def produce_loaded_user_covmat_path( @configparser.explicit_node def produce_nnfit_theory_covmat( - self, - use_thcovmat_in_sampling: bool, - use_thcovmat_in_fitting: bool, - inclusive_use_scalevar_uncertainties, - use_user_uncertainties: bool = False, + self, inclusive_use_scalevar_uncertainties, use_user_uncertainties: bool = False ): """ Return the theory covariance matrix used in the fit. + + This function is only used in vp-setupfit to store the necessary covmats as .csv files in + the tables directory. """ if inclusive_use_scalevar_uncertainties: if use_user_uncertainties: @@ -1153,13 +1152,7 @@ def produce_nnfit_theory_covmat( f = user_covmat_fitting - @functools.wraps(f) - def res(*args, **kwargs): - return f(*args, **kwargs) - - # Set this to get the same filename regardless of the action. - res.__name__ = "theory_covmat" - return res + return f def produce_fitthcovmat( self, use_thcovmat_if_present: bool = False, fit: (str, type(None)) = None @@ -1630,21 +1623,6 @@ def produce_group_dataset_inputs_by_metadata(self, data_input, processed_metadat for name, group in res.items() ] - def produce_fivetheories(self, point_prescription): - if point_prescription == "5bar point": - return "bar" - elif point_prescription == "5 point": - return "nobar" - return None - - def produce_seventheories(self, point_prescription): - if point_prescription == "7 point": - # This is None because is the default choice - return None - elif point_prescription == "7original point": - return "original" - return None - def produce_group_dataset_inputs_by_experiment(self, data_input): return self.produce_group_dataset_inputs_by_metadata(data_input, "experiment") diff --git a/validphys2/src/validphys/theorycovariance/construction.py b/validphys2/src/validphys/theorycovariance/construction.py index 4e9a200f22..ea8aab30ef 100644 --- a/validphys2/src/validphys/theorycovariance/construction.py +++ b/validphys2/src/validphys/theorycovariance/construction.py @@ -27,12 +27,7 @@ @check_using_theory_covmat def theory_covmat_dataset( - results, - results_central_bytheoryids, - use_theorycovmat, # for the check - point_prescription, - fivetheories=None, - seventheories=None, + results, results_central_bytheoryids, point_prescription, use_theorycovmat # for the check ): """ Compute the theory covmat for a collection of theoryids for a single dataset. @@ -50,9 +45,7 @@ def theory_covmat_dataset( # Compute the theory contribution to the covmats deltas = list(t.central_value - cv for t in theory_results) - thcovmat = compute_covs_pt_prescrip( - point_prescription, l, "A", deltas, fivetheories=fivetheories, seventheories=seventheories - ) + thcovmat = compute_covs_pt_prescrip(point_prescription, l, "A", deltas) return thcovmat @@ -242,16 +235,7 @@ def covmat_n3lo_ad(name1, name2, deltas1, deltas2): return 1 / norm * s -def compute_covs_pt_prescrip( - point_prescription, - l, - name1, - deltas1, - name2=None, - deltas2=None, - fivetheories=None, - seventheories=None, -): +def compute_covs_pt_prescrip(point_prescription, l, name1, deltas1, name2=None, deltas2=None): """Utility to compute the covariance matrix by prescription given the shifts with respect to the central value for a pair of processes. @@ -275,10 +259,6 @@ def compute_covs_pt_prescrip( Process name of the second set of shifts deltas2: list(np.ndarray) list of shifts for each of the non-central theories - fivetheories: str - 5-point prescription variation - seventheories: str - 7-point prescription variation """ if name2 is None and deltas2 is not None: raise ValueError( @@ -298,22 +278,19 @@ def compute_covs_pt_prescrip( s = covmat_3fpt(name1, name2, deltas1, deltas2) elif point_prescription == "3r point": s = covmat_3rpt(name1, name2, deltas1, deltas2) - else: + elif point_prescription == "3 point": s = covmat_3pt(name1, name2, deltas1, deltas2) elif l == 5: - # 5 point -------------------------------------------------------------- - if fivetheories == "nobar": + if point_prescription == "5 point": s = covmat_5pt(name1, name2, deltas1, deltas2) - # 5bar-point ----------------------------------------------------------- - else: + elif point_prescription == "5bar point": s = covmat_5barpt(name1, name2, deltas1, deltas2) - # --------------------------------------------------------------------- elif l == 7: # Outdated 7pts implementation: left for posterity --------------------- - if seventheories == "original": + if point_prescription == "7original point": s = covmat_7pt_orig(name1, name2, deltas1, deltas2) # 7pt (Gavin) ---------------------------------------------------------- - else: + elif point_prescription == "7 point": s = covmat_7pt(name1, name2, deltas1, deltas2) elif l == 9: s = covmat_9pt(name1, name2, deltas1, deltas2) @@ -361,7 +338,7 @@ def compute_covs_pt_prescrip( @check_correct_theory_combination -def covs_pt_prescrip(combine_by_type, theoryids, point_prescription, fivetheories, seventheories): +def covs_pt_prescrip(combine_by_type, theoryids, point_prescription): """Produces the sub-matrices of the theory covariance matrix according to a point prescription which matches the number of input theories. If 5 theories are provided, a scheme 'bar' or 'nobar' must be @@ -386,14 +363,7 @@ def covs_pt_prescrip(combine_by_type, theoryids, point_prescription, fivetheorie central2, *others2 = process_info.preds[name2] deltas2 = list(other - central2 for other in others2) s = compute_covs_pt_prescrip( - point_prescription, - len(theoryids), - name1, - deltas1, - name2, - deltas2, - fivetheories, - seventheories, + point_prescription, len(theoryids), name1, deltas1, name2, deltas2 ) start_locs = (start_proc[name1], start_proc[name2]) covmats[start_locs] = s diff --git a/validphys2/src/validphys/theorycovariance/output.py b/validphys2/src/validphys/theorycovariance/output.py index 6df269ad2d..3501d0cf54 100644 --- a/validphys2/src/validphys/theorycovariance/output.py +++ b/validphys2/src/validphys/theorycovariance/output.py @@ -2,7 +2,6 @@ output.py Basic tools for plotting theory covariance matrices and their properties. """ -from __future__ import generator_stop import logging from math import inf @@ -195,22 +194,22 @@ def plot_expcorrmat_heatmap(procs_corrmat): @figure -def plot_normthcovmat_heatmap_custom(theory_normcovmat_custom, theoryids, fivetheories): +def plot_normthcovmat_heatmap_custom(theory_normcovmat_custom, theoryids, point_prescription): """Matrix plot for block diagonal theory covariance matrix by process type""" l = len(theoryids) if l == 5: - if fivetheories == "bar": + if point_prescription == "5bar point": l = r"$\bar{5}$" fig = plot_covmat_heatmap(theory_normcovmat_custom, f"Theory Covariance matrix ({l} pt)") return fig @figure -def plot_thcorrmat_heatmap_custom(theory_corrmat_custom, theoryids, fivetheories): +def plot_thcorrmat_heatmap_custom(theory_corrmat_custom, theoryids, point_prescription): """Matrix plot of the theory correlation matrix, correlations by process type""" l = len(theoryids) if l == 5: - if fivetheories == "bar": + if point_prescription == "5bar point": l = r"$\bar{5}$" fig = plot_corrmat_heatmap(theory_corrmat_custom, f"Theory Correlation matrix ({l} pt)") return fig @@ -218,12 +217,12 @@ def plot_thcorrmat_heatmap_custom(theory_corrmat_custom, theoryids, fivetheories @figure def plot_expplusthcorrmat_heatmap_custom( - experimentplustheory_corrmat_custom, theoryids, fivetheories + experimentplustheory_corrmat_custom, theoryids, point_prescription ): """Matrix plot of the exp + theory correlation matrix""" l = len(theoryids) if l == 5: - if fivetheories == "bar": + if point_prescription == "5bar point": l = r"$\bar{5}$" fig = plot_corrmat_heatmap( experimentplustheory_corrmat_custom, f"Experimental + Theory Correlation Matrix ({l} pt)" @@ -233,12 +232,12 @@ def plot_expplusthcorrmat_heatmap_custom( @figure def plot_diag_cov_comparison( - theory_covmat_custom, procs_covmat, procs_data_values, theoryids, fivetheories + theory_covmat_custom, procs_covmat, procs_data_values, theoryids, point_prescription ): """Plot of sqrt(cov_ii)/|data_i| for cov = exp, theory, exp+theory""" l = len(theoryids) if l == 5: - if fivetheories == "bar": + if point_prescription == "5bar point": l = r"$\bar{5}$" data = np.abs(procs_data_values) sqrtdiags_th = np.sqrt(np.diag(theory_covmat_custom)) / data diff --git a/validphys2/src/validphys/theorycovariance/tests.py b/validphys2/src/validphys/theorycovariance/tests.py index ae6d1112ce..53de9083ee 100644 --- a/validphys2/src/validphys/theorycovariance/tests.py +++ b/validphys2/src/validphys/theorycovariance/tests.py @@ -1,9 +1,7 @@ -# -*- coding: utf-8 -*- """ tests.py Tools for testing theory covariance matrices and their properties. """ -from __future__ import generator_stop from collections import namedtuple import logging @@ -36,9 +34,7 @@ def shift_vector(matched_dataspecs_results, process, dataset_name): norm_shifts = shifts / norm dsnames = np.full(len(shifts), dataset_name, dtype=object) processnames = np.full(len(shifts), process, dtype=object) - index = pd.MultiIndex.from_arrays( - [processnames, dsnames], names=["group", "dataset"] - ) + index = pd.MultiIndex.from_arrays([processnames, dsnames], names=["group", "dataset"]) return pd.DataFrame({'shifts': norm_shifts, 'norm': norm}, index=index) @@ -296,9 +292,8 @@ def ordered_alltheory_vector(alltheory_vector, doubleindex_thcovmat): def evals_nonzero_basis( ordered_alltheory_vector, theory_covmat_custom, - fivetheories, + point_prescription, orthonormalisation: (str, type(None)) = None, - seventheories: (str, type(None)) = None, ): covmat = theory_covmat_custom / ( np.outer(ordered_alltheory_vector[0], ordered_alltheory_vector[0]) @@ -343,11 +338,11 @@ def evals_nonzero_basis( # (mu_0; mu_1, mu_2, ..., mu_p) if num_pts == 3: xs = vectors_3pt(splitdiffs) - elif (num_pts == 5) and (fivetheories == "nobar"): + elif (num_pts == 5) and (point_prescription == "5 point"): xs = vectors_5pt(splitdiffs) - elif (num_pts == 5) and (fivetheories == "bar"): + elif (num_pts == 5) and (point_prescription == "5bar point"): xs = vectors_5barpt(splitdiffs) - elif (num_pts == 7) and (seventheories != "original"): + elif (num_pts == 7) and (point_prescription != "7original point"): xs = vectors_7pt(splitdiffs) elif num_pts == 9: xs = vectors_9pt(splitdiffs) @@ -381,11 +376,13 @@ def projected_condition_num(evals_nonzero_basis): cond_num = evals_nonzero_basis[2] return cond_num + def fnorm_shifts_ordered(concatenated_shx_vector, doubleindex_thcovmat): """Shifts vectors ordered as the theory covmat.""" fnorm_vector = fnorm_shifts_byprocess(concatenated_shx_vector, doubleindex_thcovmat) return np.array(sum(fnorm_vector, [])) + def theory_shift_test(fnorm_shifts_ordered, evals_nonzero_basis): """Compares the NNLO-NLO shift, f, with the eigenvectors and eigenvalues of the theory covariance matrix, and returns the component of the NNLO-NLO shift @@ -501,7 +498,9 @@ def projector_eigenvalue_ratio(theory_shift_test): @figure -def eigenvector_plot(evals_nonzero_basis, fnorm_shifts_ordered, tripleindex_thcovmat_complete, dataspecs): +def eigenvector_plot( + evals_nonzero_basis, fnorm_shifts_ordered, tripleindex_thcovmat_complete, dataspecs +): """Produces a plot of the eigenvectors for the projected matrix, transformed back to the data space.""" evals = evals_nonzero_basis[0][::-1] @@ -518,7 +517,11 @@ def eigenvector_plot(evals_nonzero_basis, fnorm_shifts_ordered, tripleindex_thco eval_3sf = floatformatting.significant_digits(eval.item(), 3) evec = pd.DataFrame(evec, index=tripleindex_thcovmat_complete) evec = evec.reindex(newindex) - ax.plot(-1.0 * (f.values), color="k", label=f"{dataspecs[1]['speclabel']}-{dataspecs[0]['speclabel']} shift") + ax.plot( + -1.0 * (f.values), + color="k", + label=f"{dataspecs[1]['speclabel']}-{dataspecs[0]['speclabel']} shift", + ) ax.plot(evec.values, label="Eigenvector") ticklocs, ticklabels, startlocs = matrix_plot_labels(evec) # Shift startlocs elements 0.5 to left so lines are between indexes @@ -558,7 +561,12 @@ def deltamiss_plot( fmiss_reordered = fmiss.reindex(f.index) # Plotting fig, ax = plotutils.subplots(figsize=(20, 10)) - ax.plot(f.values * 100, ".-", label=f"{dataspecs[1]['speclabel']}-{dataspecs[0]['speclabel']} Shift", color="black") + ax.plot( + f.values * 100, + ".-", + label=f"{dataspecs[1]['speclabel']}-{dataspecs[0]['speclabel']} Shift", + color="black", + ) ax.plot( fmiss_reordered.values * 100, ".-", label=r"$\delta_{miss}$" + f" ({l} pt)", color="blue" ) @@ -579,7 +587,9 @@ def deltamiss_plot( def diagdf_theory_covmat(theory_covmat_custom): """Return a Dataframe indexed with groups and dataset of the diagonal entry of the theory covmat.""" - return pd.DataFrame(data=np.diag(theory_covmat_custom.values), index=theory_covmat_custom.index.droplevel(2)) + return pd.DataFrame( + data=np.diag(theory_covmat_custom.values), index=theory_covmat_custom.index.droplevel(2) + ) def doubleindex_set_byprocess(group_dataset_inputs_by_process): @@ -596,17 +606,19 @@ def concatenated_shx_vector(shx_vector): return pd.concat(shx_vector) -def sqrtdiags_thcovmat_byprocess(doubleindex_set_byprocess, diagdf_theory_covmat, concatenated_shx_vector): +def sqrtdiags_thcovmat_byprocess( + doubleindex_set_byprocess, diagdf_theory_covmat, concatenated_shx_vector +): """Ratio of the sqrts of the diagonal entries of the theory covmat to the normalization entries of the shift vectors, ordered by process.""" sqrtdiags = [] for index in doubleindex_set_byprocess: if isinstance(concatenated_shx_vector.loc[index[0]].loc[index[1]].norm, np.float64): sqrtdiags.append( - list( - np.sqrt(diagdf_theory_covmat.loc[index[0]].loc[index[1]].values) - / concatenated_shx_vector.loc[index[0]].loc[index[1]].norm + list( + np.sqrt(diagdf_theory_covmat.loc[index[0]].loc[index[1]].values) + / concatenated_shx_vector.loc[index[0]].loc[index[1]].norm + ) ) - ) else: sqrtdiags.append( list( @@ -634,16 +646,18 @@ def ticklocs_thcovmat(theory_covmat_custom): @figure def shift_diag_cov_comparison( - sqrtdiags_thcovmat_byprocess, fnorm_shifts_byprocess, point_prescription, ticklocs_thcovmat, dataspecs + sqrtdiags_thcovmat_byprocess, + fnorm_shifts_byprocess, + point_prescription, + ticklocs_thcovmat, + dataspecs, ): """Plot of the comparison of a shift between two pertubative order and the diagonal entries of the theory covmat, both normalized to the first of the two perturbative orders.""" # Concatenation of the arrays fnorm_concat = sum(fnorm_shifts_byprocess, []) sqrtdiags_concat = sum(sqrtdiags_thcovmat_byprocess, []) fig, ax = plotutils.subplots(figsize=(20, 10)) - ax.plot( - np.array(sqrtdiags_concat), ".-", label=f"MHOU ({point_prescription})", color="red" - ) + ax.plot(np.array(sqrtdiags_concat), ".-", label=f"MHOU ({point_prescription})", color="red") ax.plot(-np.array(sqrtdiags_concat), ".-", color="red") ax.plot( -np.array(fnorm_concat), diff --git a/validphys2/src/validphys/theorycovariance/theorycovarianceutils.py b/validphys2/src/validphys/theorycovariance/theorycovarianceutils.py index 7ff8628c23..11ec428f9a 100644 --- a/validphys2/src/validphys/theorycovariance/theorycovarianceutils.py +++ b/validphys2/src/validphys/theorycovariance/theorycovarianceutils.py @@ -14,7 +14,7 @@ def check_correct_theory_combination_internal( - theoryids, fivetheories, point_prescription: (str, type(None)) = None + theoryids, point_prescription: (str, type(None)) = None ): """Checks that a valid theory combination corresponding to an existing prescription has been inputted""" @@ -40,18 +40,10 @@ def check_correct_theory_combination_internal( correct_xifs = [1.0, 2.0, 0.5] correct_xirs = [1.0, 2.0, 0.5] elif l == 5: - check( - fivetheories is not None, - "For five input theories a prescription bar or nobar" - "for the flag fivetheories must be specified.", - ) - check( - fivetheories in opts, "Invalid choice of prescription for 5 points", fivetheories, opts - ) - if fivetheories == "nobar": + if point_prescription == "5 point": correct_xifs = [1.0, 2.0, 0.5, 1.0, 1.0] correct_xirs = [1.0, 1.0, 1.0, 2.0, 0.5] - elif fivetheories == "bar": + elif point_prescription == "5bar point": correct_xifs = [1.0, 2.0, 0.5, 2.0, 0.5] correct_xirs = [1.0, 2.0, 0.5, 0.5, 2.0] elif l == 7: