comparison COBRAxy/ras_to_bounds_beta.py @ 414:5086145cfb96 draft

Uploaded
author francesco_lapi
date Mon, 08 Sep 2025 21:54:14 +0000
parents 6b015d3184ab
children 5f8f4a2d1370
comparison
equal deleted inserted replaced
413:7a3ccf066b2c 414:5086145cfb96
8 from cobra import Model, Reaction, Metabolite 8 from cobra import Model, Reaction, Metabolite
9 import re 9 import re
10 import sys 10 import sys
11 import csv 11 import csv
12 from joblib import Parallel, delayed, cpu_count 12 from joblib import Parallel, delayed, cpu_count
13 import utils.rule_parsing as rulesUtils
13 14
14 # , medium 15 # , medium
15 16
16 ################################# process args ############################### 17 ################################# process args ###############################
17 def process_args(args :List[str] = None) -> argparse.Namespace: 18 def process_args(args :List[str] = None) -> argparse.Namespace:
147 if upper_bound!=0 and lower_bound!=0: 148 if upper_bound!=0 and lower_bound!=0:
148 new_bounds.loc[reaction, "lower_bound"] = valMin 149 new_bounds.loc[reaction, "lower_bound"] = valMin
149 new_bounds.loc[reaction, "upper_bound"] = valMax 150 new_bounds.loc[reaction, "upper_bound"] = valMax
150 return new_bounds 151 return new_bounds
151 152
153 ################################- DATA GENERATION -################################
154 ReactionId = str
155 def generate_rules(model: cobra.Model, *, asParsed = True) -> Union[Dict[ReactionId, rulesUtils.OpList], Dict[ReactionId, str]]:
156 """
157 Generates a dictionary mapping reaction ids to rules from the model.
158
159 Args:
160 model : the model to derive data from.
161 asParsed : if True parses the rules to an optimized runtime format, otherwise leaves them as strings.
162
163 Returns:
164 Dict[ReactionId, rulesUtils.OpList] : the generated dictionary of parsed rules.
165 Dict[ReactionId, str] : the generated dictionary of raw rules.
166 """
167 # Is the below approach convoluted? yes
168 # Ok but is it inefficient? probably
169 # Ok but at least I don't have to repeat the check at every rule (I'm clinically insane)
170 _ruleGetter = lambda reaction : reaction.gene_reaction_rule
171 ruleExtractor = (lambda reaction :
172 rulesUtils.parseRuleToNestedList(_ruleGetter(reaction))) if asParsed else _ruleGetter
173
174 return {
175 reaction.id : ruleExtractor(reaction)
176 for reaction in model.reactions
177 if reaction.gene_reaction_rule }
178
179 def generate_reactions(model :cobra.Model, *, asParsed = True) -> Dict[ReactionId, str]:
180 """
181 Generates a dictionary mapping reaction ids to reaction formulas from the model.
182
183 Args:
184 model : the model to derive data from.
185 asParsed : if True parses the reactions to an optimized runtime format, otherwise leaves them as they are.
186
187 Returns:
188 Dict[ReactionId, str] : the generated dictionary.
189 """
190
191 unparsedReactions = {
192 reaction.id : reaction.reaction
193 for reaction in model.reactions
194 if reaction.reaction
195 }
196
197 if not asParsed: return unparsedReactions
198
199 return reactionUtils.create_reaction_dict(unparsedReactions)
200
201 def get_medium(model:cobra.Model) -> pd.DataFrame:
202 trueMedium=[]
203 for r in model.reactions:
204 positiveCoeff=0
205 for m in r.metabolites:
206 if r.get_coefficient(m.id)>0:
207 positiveCoeff=1;
208 if (positiveCoeff==0 and r.lower_bound<0):
209 trueMedium.append(r.id)
210
211 df_medium = pd.DataFrame()
212 df_medium["reaction"] = trueMedium
213 return df_medium
214
215 def generate_bounds(model:cobra.Model) -> pd.DataFrame:
216
217 rxns = []
218 for reaction in model.reactions:
219 rxns.append(reaction.id)
220
221 bounds = pd.DataFrame(columns = ["lower_bound", "upper_bound"], index=rxns)
222
223 for reaction in model.reactions:
224 bounds.loc[reaction.id] = [reaction.lower_bound, reaction.upper_bound]
225 return bounds
226
227
228
229 def generate_compartments(model: cobra.Model) -> pd.DataFrame:
230 """
231 Generates a DataFrame containing compartment information for each reaction.
232 Creates columns for each compartment position (Compartment_1, Compartment_2, etc.)
233
234 Args:
235 model: the COBRA model to extract compartment data from.
236
237 Returns:
238 pd.DataFrame: DataFrame with ReactionID and compartment columns
239 """
240 pathway_data = []
241
242 # First pass: determine the maximum number of pathways any reaction has
243 max_pathways = 0
244 reaction_pathways = {}
245
246 for reaction in model.reactions:
247 # Get unique pathways from all metabolites in the reaction
248 if type(reaction.annotation['pathways']) == list:
249 reaction_pathways[reaction.id] = reaction.annotation['pathways']
250 max_pathways = max(max_pathways, len(reaction.annotation['pathways']))
251 else:
252 reaction_pathways[reaction.id] = [reaction.annotation['pathways']]
253
254 # Create column names for pathways
255 pathway_columns = [f"Pathway_{i+1}" for i in range(max_pathways)]
256
257 # Second pass: create the data
258 for reaction_id, pathways in reaction_pathways.items():
259 row = {"ReactionID": reaction_id}
260
261 # Fill pathway columns
262 for i in range(max_pathways):
263 col_name = pathway_columns[i]
264 if i < len(pathways):
265 row[col_name] = pathways[i]
266 else:
267 row[col_name] = None # or "" if you prefer empty strings
268
269 pathway_data.append(row)
270
271 return pd.DataFrame(pathway_data)
272
152 def save_model(model, filename, output_folder, file_format='csv'): 273 def save_model(model, filename, output_folder, file_format='csv'):
153 """ 274 """
154 Save a COBRA model to file in the specified format. 275 Save a COBRA model to file in the specified format.
155 276
156 Args: 277 Args:
168 try: 289 try:
169 if file_format == 'tabular' or file_format == 'csv': 290 if file_format == 'tabular' or file_format == 'csv':
170 # Special handling for tabular format using utils functions 291 # Special handling for tabular format using utils functions
171 filepath = os.path.join(output_folder, f"{filename}.csv") 292 filepath = os.path.join(output_folder, f"{filename}.csv")
172 293
173 rules = utils.generate_rules(model, asParsed = False) 294 rules = generate_rules(model, asParsed = False)
174 reactions = utils.generate_reactions(model, asParsed = False) 295 reactions = generate_reactions(model, asParsed = False)
175 bounds = utils.generate_bounds(model) 296 bounds = generate_bounds(model)
176 medium = utils.get_medium(model) 297 medium = get_medium(model)
177 298
178 try: 299 try:
179 compartments = utils.generate_compartments(model) 300 compartments = utils.generate_compartments(model)
180 except: 301 except:
181 compartments = None 302 compartments = None
267 modified_model = apply_bounds_to_model(model, new_bounds) 388 modified_model = apply_bounds_to_model(model, new_bounds)
268 save_model(modified_model, cellName, save_models_path, save_models_format) 389 save_model(modified_model, cellName, save_models_path, save_models_format)
269 390
270 pass 391 pass
271 392
272 def generate_bounds(model: cobra.Model, ras=None, output_folder='output/', save_models=False, save_models_path='saved_models/', save_models_format='csv') -> pd.DataFrame: 393 def generate_bounds_model(model: cobra.Model, ras=None, output_folder='output/', save_models=False, save_models_path='saved_models/', save_models_format='csv') -> pd.DataFrame:
273 """ 394 """
274 Generate reaction bounds for a metabolic model based on medium conditions and optional RAS adjustments. 395 Generate reaction bounds for a metabolic model based on medium conditions and optional RAS adjustments.
275 396
276 Args: 397 Args:
277 model (cobra.Model): The metabolic model for which bounds will be generated. 398 model (cobra.Model): The metabolic model for which bounds will be generated.
367 print("\n=== VALIDAZIONE MODELLO ===") 488 print("\n=== VALIDAZIONE MODELLO ===")
368 for key, value in validation.items(): 489 for key, value in validation.items():
369 print(f"{key}: {value}") 490 print(f"{key}: {value}")
370 491
371 if(ARGS.ras_selector == True): 492 if(ARGS.ras_selector == True):
372 generate_bounds(model, ras=ras_combined, output_folder=ARGS.output_path, 493 generate_bounds_model(model, ras=ras_combined, output_folder=ARGS.output_path,
373 save_models=ARGS.save_models, save_models_path=ARGS.save_models_path, 494 save_models=ARGS.save_models, save_models_path=ARGS.save_models_path,
374 save_models_format=ARGS.save_models_format) 495 save_models_format=ARGS.save_models_format)
375 class_assignments.to_csv(ARGS.cell_class, sep='\t', index=False) 496 class_assignments.to_csv(ARGS.cell_class, sep='\t', index=False)
376 else: 497 else:
377 generate_bounds(model, output_folder=ARGS.output_path, 498 generate_bounds_model(model, output_folder=ARGS.output_path,
378 save_models=ARGS.save_models, save_models_path=ARGS.save_models_path, 499 save_models=ARGS.save_models, save_models_path=ARGS.save_models_path,
379 save_models_format=ARGS.save_models_format) 500 save_models_format=ARGS.save_models_format)
380 501
381 pass 502 pass
382 503