Mercurial > repos > bimib > cobraxy
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 |
