Mercurial > repos > bimib > cobraxy

comparison COBRAxy/utils/ras_to_bounds.py @ 121:2121cc1aef45 draft

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Uploaded
author luca_milaz
date Mon, 14 Oct 2024 08:19:38 +0000
parents
children
comparison
equal deleted inserted replaced
120:567ffd4333d4 121:2121cc1aef45
1 import argparse
2 import utils.general_utils as utils
3 from typing import Optional, List
4 import os
5 import numpy as np
6 import pandas as pd
7 import cobra
8 import sys
9 import csv
10 from joblib import Parallel, delayed, cpu_count
11
12 ################################# process args ###############################
13 def process_args(args :List[str]) -> argparse.Namespace:
14 """
15 Processes command-line arguments.
16
17 Args:
18 args (list): List of command-line arguments.
19
20 Returns:
21 Namespace: An object containing parsed arguments.
22 """
23 parser = argparse.ArgumentParser(usage = '%(prog)s [options]',
24 description = 'process some value\'s')
25
26 parser.add_argument(
27 '-ms', '--model_selector',
28 type = utils.Model, default = utils.Model.ENGRO2, choices = [utils.Model.ENGRO2, utils.Model.Custom],
29 help = 'chose which type of model you want use')
30
31 parser.add_argument("-mo", "--model", type = str,
32 help = "path to input file with custom rules, if provided")
33
34 parser.add_argument("-mn", "--model_name", type = str, help = "custom mode name")
35
36 parser.add_argument(
37 '-mes', '--medium_selector',
38 default = "allOpen",
39 help = 'chose which type of medium you want use')
40
41 parser.add_argument("-meo", "--medium", type = str,
42 help = "path to input file with custom medium, if provided")
43
44 parser.add_argument('-ol', '--out_log',
45 help = "Output log")
46
47 parser.add_argument('-td', '--tool_dir',
48 type = str,
49 required = True,
50 help = 'your tool directory')
51
52 parser.add_argument('-ir', '--input_ras',
53 type=str,
54 required = False,
55 help = 'input ras')
56
57 parser.add_argument('-rn', '--name',
58 type=str,
59 help = 'ras class names')
60
61 parser.add_argument('-rs', '--ras_selector',
62 required = True,
63 type=utils.Bool("using_RAS"),
64 help = 'ras selector')
65
66 parser.add_argument('-cc', '--cell_class',
67 type = str,
68 help = 'output of cell class')
69
70
71 ARGS = parser.parse_args()
72 return ARGS
73
74 ########################### warning ###########################################
75 def warning(s :str) -> None:
76 """
77 Log a warning message to an output log file and print it to the console.
78
79 Args:
80 s (str): The warning message to be logged and printed.
81
82 Returns:
83 None
84 """
85 with open(ARGS.out_log, 'a') as log:
86 log.write(s + "\n\n")
87 print(s)
88
89 ############################ dataset input ####################################
90 def read_dataset(data :str, name :str) -> pd.DataFrame:
91 """
92 Read a dataset from a CSV file and return it as a pandas DataFrame.
93
94 Args:
95 data (str): Path to the CSV file containing the dataset.
96 name (str): Name of the dataset, used in error messages.
97
98 Returns:
99 pandas.DataFrame: DataFrame containing the dataset.
100
101 Raises:
102 pd.errors.EmptyDataError: If the CSV file is empty.
103 sys.exit: If the CSV file has the wrong format, the execution is aborted.
104 """
105 try:
106 dataset = pd.read_csv(data, sep = '\t', header = 0, engine='python')
107 except pd.errors.EmptyDataError:
108 sys.exit('Execution aborted: wrong format of ' + name + '\n')
109 if len(dataset.columns) < 2:
110 sys.exit('Execution aborted: wrong format of ' + name + '\n')
111 return dataset
112
113
114 def apply_ras_bounds(model, ras_row, rxns_ids, mediumRxns_ids):
115 """
116 Adjust the bounds of reactions in the model based on RAS values.
117
118 Args:
119 model (cobra.Model): The metabolic model to be modified.
120 ras_row (pd.Series): A row from a RAS DataFrame containing scaling factors for reaction bounds.
121 rxns_ids (list of str): List of reaction IDs to which the scaling factors will be applied.
122 mediumRxns_ids (list of str): List of reaction IDs in the medium. Their RAS is set to zero, but they are already set in the model.
123 Returns:
124 None
125 """
126 for reaction in ras_row.index:
127 scaling_factor = ras_row[reaction]
128 if(scaling_factor not in [np.nan, None]):
129 lower_bound=model.reactions.get_by_id(reaction).lower_bound
130 upper_bound=model.reactions.get_by_id(reaction).upper_bound
131 valMax=float((upper_bound)*scaling_factor)
132 valMin=float((lower_bound)*scaling_factor)
133 if upper_bound!=0 and lower_bound==0:
134 model.reactions.get_by_id(reaction).upper_bound=valMax
135 if upper_bound==0 and lower_bound!=0:
136 model.reactions.get_by_id(reaction).lower_bound=valMin
137 if upper_bound!=0 and lower_bound!=0:
138 model.reactions.get_by_id(reaction).lower_bound=valMin
139 model.reactions.get_by_id(reaction).upper_bound=valMax
140 pass
141
142 def process_ras_cell(cellName, ras_row, model, rxns_ids, mediumRxns_ids, output_folder):
143 """
144 Process a single RAS cell, apply bounds, and save the bounds to a CSV file.
145
146 Args:
147 cellName (str): The name of the RAS cell (used for naming the output file).
148 ras_row (pd.Series): A row from a RAS DataFrame containing scaling factors for reaction bounds.
149 model (cobra.Model): The metabolic model to be modified.
150 rxns_ids (list of str): List of reaction IDs to which the scaling factors will be applied.
151 mediumRxns_ids (list of str): List of reaction IDs in the medium. Their RAS is set to zero, but they are already set in the model.
152 output_folder (str): Folder path where the output CSV file will be saved.
153
154 Returns:
155 None
156 """
157 model_new = model.copy()
158 apply_ras_bounds(model_new, ras_row, rxns_ids, mediumRxns_ids)
159 bounds = pd.DataFrame([(rxn.lower_bound, rxn.upper_bound) for rxn in model_new.reactions], index=rxns_ids, columns=["lower_bound", "upper_bound"])
160 bounds.to_csv(output_folder + cellName + ".csv", sep='\t', index=True)
161 pass
162
163 def generate_bounds(model: cobra.Model, medium: dict, ras=None, output_folder='output/') -> pd.DataFrame:
164 """
165 Generate reaction bounds for a metabolic model based on medium conditions and optional RAS adjustments.
166
167 Args:
168 model (cobra.Model): The metabolic model for which bounds will be generated.
169 medium (dict): A dictionary where keys are reaction IDs and values are the medium conditions.
170 ras (pd.DataFrame, optional): RAS pandas dataframe. Defaults to None.
171 output_folder (str, optional): Folder path where output CSV files will be saved. Defaults to 'output/'.
172
173 Returns:
174 pd.DataFrame: DataFrame containing the bounds of reactions in the model.
175 """
176 rxns_ids = [rxn.id for rxn in model.reactions]
177
178 # Set medium conditions
179 for reaction, value in medium.items():
180 if value is not None:
181 ## SOLO ENGRO2
182 if(reaction != "EX_thbpt_e" and reaction != "EX_lac__L_e"):
183 model.reactions.get_by_id(reaction).lower_bound = -float(value)
184 if(reaction == "EX_lac__L_e"):
185 model.reactions.get_by_id(reaction).lower_bound = float(0.0)
186
187 mediumRxns_ids = medium.keys()
188
189 # Perform Flux Variability Analysis (FVA) on this medium
190 df_FVA = cobra.flux_analysis.flux_variability_analysis(model, fraction_of_optimum=0, processes=1).round(8)
191
192 # Set FVA bounds
193 for reaction in rxns_ids:
194 model.reactions.get_by_id(reaction).lower_bound = float(df_FVA.loc[reaction, "minimum"])
195 model.reactions.get_by_id(reaction).upper_bound = float(df_FVA.loc[reaction, "maximum"])
196
197 if ras is not None:
198 #Parallel(n_jobs=cpu_count())(delayed(process_ras_cell)(cellName, ras_row, model, rxns_ids, output_folder) for cellName, ras_row in ras.iterrows())
199 for cellName, ras_row in ras.iterrows():
200 process_ras_cell(cellName, ras_row, model, rxns_ids, mediumRxns_ids, output_folder)
201 break #just one cell for testing
202 else:
203 model_new = model.copy()
204 apply_ras_bounds(model_new, pd.Series([1]*len(rxns_ids), index=rxns_ids), rxns_ids, mediumRxns_ids)
205 bounds = pd.DataFrame([(rxn.lower_bound, rxn.upper_bound) for rxn in model_new.reactions], index=rxns_ids, columns=["lower_bound", "upper_bound"])
206 bounds.to_csv(output_folder + "bounds.csv", sep='\t', index=True)
207 pass
208
209
210
211 ############################# main ###########################################
212 def main() -> None:
213 """
214 Initializes everything and sets the program in motion based on the fronted input arguments.
215
216 Returns:
217 None
218 """
219 if not os.path.exists('ras_to_bounds'):
220 os.makedirs('ras_to_bounds')
221
222
223 global ARGS
224 ARGS = process_args(sys.argv)
225
226 ARGS.output_folder = 'ras_to_bounds/'
227
228 if(ARGS.ras_selector == True):
229 ras_file_list = ARGS.input_ras.split(",")
230 ras_file_names = ARGS.name.split(",")
231 ras_class_names = []
232 for file in ras_file_names:
233 ras_class_names.append(file.split(".")[0])
234 ras_list = []
235 class_assignments = pd.DataFrame(columns=["Patient_ID", "Class"])
236 for ras_matrix, ras_class_name in zip(ras_file_list, ras_class_names):
237 ras = read_dataset(ras_matrix, "ras dataset")
238 ras.replace("None", None, inplace=True)
239 ras.set_index("Reactions", drop=True, inplace=True)
240 ras = ras.T
241 ras = ras.astype(float)
242 ras_list.append(ras)
243 for patient_id in ras.index:
244 class_assignments.loc[class_assignments.shape[0]] = [patient_id, ras_class_name]
245
246
247 # Concatenate all ras DataFrames into a single DataFrame
248 ras_combined = pd.concat(ras_list, axis=0)
249 # Normalize the RAS values by max RAS
250 ras_combined = ras_combined.div(ras_combined.max(axis=0))
251 #ras_combined = ras_combined.fillna(0)
252 #il ras c'è per tutti o non c'è per nessuno
253
254
255
256 model_type :utils.Model = ARGS.model_selector
257 if model_type is utils.Model.Custom:
258 model = model_type.getCOBRAmodel(customPath = utils.FilePath.fromStrPath(ARGS.model), customExtension = utils.FilePath.fromStrPath(ARGS.model_name).ext)
259 else:
260 model = model_type.getCOBRAmodel(toolDir=ARGS.tool_dir)
261
262 if(ARGS.medium_selector == "Custom"):
263 medium = read_dataset(ARGS.medium, "medium dataset")
264 medium.set_index(medium.columns[0], inplace=True)
265 medium = medium.astype(float)
266 medium = medium[medium.columns[0]].to_dict()
267 else:
268 df_mediums = pd.read_csv(ARGS.tool_dir + "/local/medium/medium.csv", index_col = 0)
269 ARGS.medium_selector = ARGS.medium_selector.replace("_", " ")
270 medium = df_mediums[[ARGS.medium_selector]]
271 medium = medium[ARGS.medium_selector].to_dict()
272
273 if(ARGS.ras_selector == True):
274 generate_bounds(model, medium, ras = ras_combined, output_folder=ARGS.output_folder)
275 class_assignments.to_csv(ARGS.cell_class, sep = '\t', index = False)
276 else:
277 generate_bounds(model, medium, output_folder=ARGS.output_folder)
278
279 pass
280
281 ##############################################################################
282 if __name__ == "__main__":
283 main()