From e61ef3b8b7f7f6725d8abf5b111c2238e4a47643 Mon Sep 17 00:00:00 2001 From: SotirisTouliopoulos Date: Wed, 14 Aug 2024 00:36:39 +0300 Subject: [PATCH] correction to correlation labels output --- dingo/utils.py | 43 ++++++++++++++++++++++++++++--------------- tests/correlation.py | 13 +++++++------ 2 files changed, 35 insertions(+), 21 deletions(-) diff --git a/dingo/utils.py b/dingo/utils.py index 473721bf..711531cf 100644 --- a/dingo/utils.py +++ b/dingo/utils.py @@ -204,13 +204,14 @@ def get_matrices_of_full_dim_polytope(A, b, Aeq, beq): -def correlated_reactions(steady_states, pearson_cutoff = 0.90, indicator_cutoff = 10, +def correlated_reactions(steady_states, reactions=[], pearson_cutoff = 0.90, indicator_cutoff = 10, cells = 10, cop_coeff = 0.3, lower_triangle = True): """A Python function to calculate the pearson correlation matrix of a model and filter values based on the copula's indicator Keyword arguments: - steady_states -- A numpy array of the generated steady states fluxes + steady_states -- A numpy array of the generated steady states fluxes + reactions -- A list with the model's reactions pearson_cutoff -- A cutoff to filter reactions based on pearson coefficient indicator_cutoff -- A cutoff to filter reactions based on indicator value cells -- Number of cells to compute the copula @@ -256,6 +257,10 @@ def correlated_reactions(steady_states, pearson_cutoff = 0.90, indicator_cutoff np.fill_diagonal(filtered_corr_matrix, 1) return filtered_corr_matrix else: + # a dictionary that will store for each filtered reaction combination, + # the pearson correlation value, the copula's indicator value + # and the correlation classification + indicator_dict = {} # keep only the lower triangle corr_matrix = np.tril(corr_matrix) @@ -264,17 +269,16 @@ def correlated_reactions(steady_states, pearson_cutoff = 0.90, indicator_cutoff # find indices of correlation matrix where correlation occurs corr_indices = np.argwhere((corr_matrix > pearson_cutoff) | (corr_matrix < -pearson_cutoff)) - - # count reactions with positive or negative correlation based on indicator - positive = 0 - negative = 0 # compute copula for each set of correlated reactions for i in range(0, corr_indices.shape[0]): index1 = corr_indices[i][0] index2 = corr_indices[i][1] - + + reaction1 = reactions[index1] + reaction2 = reactions[index2] + flux1 = steady_states[index1] flux2 = steady_states[index2] @@ -302,26 +306,35 @@ def correlated_reactions(steady_states, pearson_cutoff = 0.90, indicator_cutoff # classify specific pair of reactions as positive or negative # correlated based on indicator cutoff if indicator > indicator_cutoff: - positive += 1 + pearson = filtered_corr_matrix[index1, index2] + indicator_dict[reaction1 + "_" + reaction2] = {'pearson': pearson, + 'indicator': indicator, + 'classification': "positive"} + elif indicator < 1/indicator_cutoff: - negative += 1 + pearson = filtered_corr_matrix[index1, index2] + indicator_dict[reaction1 + "_" + reaction2] = {'pearson': pearson, + 'indicator': indicator, + 'classification': "negative"} + # if they do not overcome the cutoff replace their corresponding # value in the correlation matrix with 0 else: filtered_corr_matrix[index1, index2] = 0 filtered_corr_matrix[index2, index1] = 0 + pearson = filtered_corr_matrix[index1, index2] + indicator_dict[reaction1 + "_" + reaction2] = {'pearson': pearson, + 'indicator': indicator, + 'classification': "no correlation"} print("Completed process of",i+1,"from",corr_indices.shape[0],"copulas") - - - print(positive, "out of", i+1, "copulas were positive correlated based on copula indicator") - print(negative, "out of", i+1, "copulas were negative correlated based on copula indicator") if lower_triangle == True: filtered_corr_matrix[np.triu_indices(filtered_corr_matrix.shape[0], 1)] = np.nan np.fill_diagonal(filtered_corr_matrix, 1) - return filtered_corr_matrix + return filtered_corr_matrix, indicator_dict else: np.fill_diagonal(filtered_corr_matrix, 1) - return filtered_corr_matrix \ No newline at end of file + return filtered_corr_matrix, indicator_dict + \ No newline at end of file diff --git a/tests/correlation.py b/tests/correlation.py index 88239dd0..0e935bb0 100644 --- a/tests/correlation.py +++ b/tests/correlation.py @@ -14,11 +14,12 @@ def test_correlation(self): sampler = PolytopeSampler(dingo_model) steady_states = sampler.generate_steady_states() - # calculate correlation matrix without filtering from copula indicator - corr_matrix = correlated_reactions(steady_states, - indicator_cutoff=2, - pearson_cutoff = 0.99999, - lower_triangle = False) + # calculate correlation matrix with filtering from copula indicator + corr_matrix, indicator_dict = correlated_reactions(steady_states, + reactions = reactions, + indicator_cutoff = 5, + pearson_cutoff = 0.999999, + lower_triangle = False) # sum values in the diagonal of the correlation matrix ==> 95*pearson ==> 95*1 self.assertTrue(np.trace(corr_matrix) == len(reactions)) @@ -35,7 +36,7 @@ def test_correlation(self): # calculate correlation matrix without filtering from copula indicator corr_matrix = correlated_reactions(steady_states, - indicator_cutoff=0, + indicator_cutoff = 0, pearson_cutoff = 0, lower_triangle = True)