order eigenvalues and eigenvectors

Author: andrpie

nnpdf_data/nnpdf_data/commondata/ATLAS_WJ_8TEV/filter_utils.py

@@ -5,7 +5,7 @@
 import pandas as pd
 import yaml
 
-from nnpdf_data.filter_utils.utils import matlist_to_matrix, prettify_float, symmetrize_errors
+from nnpdf_data.filter_utils.utils import matlist_to_matrix, prettify_float, symmetrize_errors, decompose_covmat
 
 yaml.add_representer(float, prettify_float)
 
@@ -406,8 +406,8 @@ def generate_data(self, variant='default', save_to_yaml=False, path='./'):
 
         # Get statistical (artidicial uncertainties)
         stat_covmat = self.__build_abs_stat_covmat()
-        eigvals, eigvecs = np.linalg.eig(stat_covmat)
-        art_stat = np.sqrt(eigvals) * eigvecs
+        
+        art_stat = decompose_covmat(stat_covmat)
 
         sys_artificial = []  # Initialize vector of artificial uncertainties
 

nnpdf_data/nnpdf_data/commondata/ATLAS_Z0_8TEV_ZMASS/filter.py

@@ -2,7 +2,7 @@
 import numpy as np
 import yaml
 
-from nnpdf_data.filter_utils.utils import prettify_float
+from nnpdf_data.filter_utils.utils import prettify_float, decompose_covmat
 
 yaml.add_representer(float, prettify_float)
 
@@ -50,8 +50,8 @@ def filter_ATLAS_Z0_8TEV_uncertainties():
 
     # compute decomposition of covariance matrix so as to get artificial systematics
     # TODO: use utils once merged in master
-    lamb, mat = np.linalg.eig(cov_matrix)
-    art_sys = np.multiply(np.sqrt(lamb), mat)
+    
+    art_sys = decompose_covmat(cov_matrix)
 
     uncertainties = []
 

nnpdf_data/nnpdf_data/commondata/CMS_WCHARM_7TEV/filter_utils.py

@@ -5,7 +5,7 @@
 import numpy as np
 import yaml
 
-from nnpdf_data.filter_utils.utils import prettify_float
+from nnpdf_data.filter_utils.utils import prettify_float, decompose_covmat
 
 yaml.add_representer(float, prettify_float)
 
@@ -211,8 +211,8 @@ def generate_data(self):
         if self.observable == 'WPWM-TOT':
             # Generate covmat and perform eigen decomposition
             covmat = self._generate_covmat(sys_unc)
-            eigvals, eigvecs = np.linalg.eig(covmat)
-            art_unc = np.sqrt(eigvals) * eigvecs
+            
+            art_unc = decompose_covmat(covmat)
 
             # Loop over bins
             for data_idx in range(len(central_data)):

nnpdf_data/nnpdf_data/commondata/CMS_WPWM_13TEV_ETA/filter_utils.py

@@ -2,7 +2,7 @@
 import uproot
 import yaml
 
-from nnpdf_data.filter_utils.utils import prettify_float
+from nnpdf_data.filter_utils.utils import prettify_float, decompose_covmat
 
 yaml.add_representer(float, prettify_float)
 
@@ -83,16 +83,6 @@ def get_data_values(version, figure):
     return data_central
 
 
-def decompose_covmat(covmat):
-    """Given a covmat it return an array sys with shape (ndat,ndat)
-    giving ndat correlated systematics for each of the ndat point.
-    The original covmat is obtained by doing sys@sys.T"""
-
-    lamb, mat = np.linalg.eig(covmat)
-    sys = np.multiply(np.sqrt(lamb), mat)
-    return sys
-
-
 def get_systematics(observable, version, figure):
     """
     Following the CMS advice we take the covariance matrix from

nnpdf_data/nnpdf_data/commondata/CMS_Z0J_8TEV/filter_utils.py

@@ -5,7 +5,7 @@
 import numpy as np
 import yaml
 
-from nnpdf_data.filter_utils.utils import prettify_float
+from nnpdf_data.filter_utils.utils import prettify_float, decompose_covmat
 
 yaml.add_representer(float, prettify_float)
 
@@ -229,8 +229,8 @@ def generate_data(self, variant='default'):
         # eigenvector basis, hence they are called "artificial uncertainties".
         # The original covmat can be reconstruted as covat = art_stat.T @ art_stat
         covmat = self._build_covmat()
-        eigvals, eigvecs = np.linalg.eig(covmat)
-        art_stat = np.sqrt(eigvals) * eigvecs * self.mult_factor
+
+        art_stat = decompose_covmat(covmat) * self.mult_factor
 
         unc_vals = []  # Initialize vector of uncertainties
         for data_idx, data in enumerate(central_data):

nnpdf_data/nnpdf_data/commondata/LHCB_Z0J_13TEV_2022/filter.py

@@ -3,7 +3,7 @@
 import numpy as np
 import yaml
 
-from nnpdf_data.filter_utils.utils import prettify_float
+from nnpdf_data.filter_utils.utils import prettify_float, decompose_covmat
 
 yaml.add_representer(float, prettify_float)
 
@@ -386,7 +386,7 @@ def processData():
                 cov = ComputeCovariance(cormat, v)
 
                 # Single value decomposition
-                sigma = OuterDecomposition(cov)
+                sigma = decompose_covmat(cov)
 
                 # Loop over the bins
                 for k in range(len(error_diag)):

nnpdf_data/nnpdf_data/commondata/LHCB_Z0_13TEV_2022/filter.py

@@ -3,7 +3,7 @@
 import numpy as np
 import yaml
 
-from nnpdf_data.filter_utils.utils import prettify_float
+from nnpdf_data.filter_utils.utils import prettify_float, decompose_covmat
 
 yaml.add_representer(float, prettify_float)
 
@@ -313,7 +313,7 @@ def processData():
                 cov = ComputeCovariance(cormat, v)
 
                 # Single value decomposition
-                sigma = OuterDecomposition(cov)
+                sigma = decompose_covmat(cov)
 
                 # Loop over the bins
                 for k in range(len(error_diag)):

nnpdf_data/nnpdf_data/filter_utils/correlations.py

@@ -1,5 +1,6 @@
 import numpy as np
 from numpy.linalg import eig
+from nnpdf_data.filter_utils.utils import sort_eigenvalues
 
 
 def upper_triangular_to_symmetric(ut, dim):
@@ -69,6 +70,7 @@ def covmat_to_artunc(ndata, covmat_list, no_of_norm_mat=0):
         b = i % ndata
         covmat[a][b] = covmat_list[i]
     eigval, eigvec = eig(covmat)
+    eigval, eigvec = sort_eigenvalues(eigval, eigvec)
     for j in range(len(eigval)):
         if eigval[j] < epsilon:
             psd_check = False

nnpdf_data/nnpdf_data/filter_utils/utils.py

@@ -106,6 +106,30 @@ def cormat_to_covmat(err_list, cormat_list):
         covmat_list.append(cormat_list[i] * err_list[a] * err_list[b])
     return covmat_list
 
+def sort_eigenvalues(evals, evecs):
+    r"""Defines ordering of the eigenvalues and eigenvectors such
+    that eigenvalues are given in decreasing order and the first
+    non-zero entry of each eigenvector is positive.
+    Parameters
+    ----------
+    evals, evecs : output of np.linalg.eig
+
+    Returns
+    -------
+    evacs_sorted, evecs_sorted : ordered eigen- values and vectors,
+    as defined above
+    
+    """
+    idx = evals.argsort()[::-1]   
+    evals_sorted = evals[idx]
+    evecs_sorted = evecs[:,idx]
+    for i in range(len(evecs_sorted)):
+        j = 0
+        while evecs_sorted[j,i] == 0:
+            j += 1
+        if evecs_sorted[j,i] < 0:
+            evecs_sorted[:, i] *= -1
+    return evals_sorted, evecs_sorted
 
 def covmat_to_artunc(ndata, covmat_list, no_of_norm_mat=0):
     r"""Convert the covariance matrix to a matrix of
@@ -152,6 +176,7 @@ def covmat_to_artunc(ndata, covmat_list, no_of_norm_mat=0):
         b = i % ndata
         covmat[a][b] = covmat_list[i]
     eigval, eigvec = eig(covmat)
+    eigval, eigvec = sort_eigenvalues(eigval, eigvec)
     for j in range(len(eigval)):
         if eigval[j] < epsilon:
             psd_check = False
@@ -397,7 +422,8 @@ def decompose_covmat(covmat):
     giving ndat correlated systematics for each of the ndat point.
     The original covmat is obtained by doing sys@sys.T"""
 
-    lamb, mat = np.linalg.eig(covmat)
+    lamb, mat = eig(covmat)
+    lamb, mat = sort_eigenvalues(lamb, mat)
     sys = np.multiply(np.sqrt(lamb), mat)
     return sys