Mercurial > repos > bimib > cobraxy

comparison COBRAxy/ras_to_bounds.py @ 127:d40d5d5435df 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:48:54 +0000
parents a9825d06cc43
children dfb6a1dc6dbc
comparison
equal deleted inserted replaced
126:30187920a779 127:d40d5d5435df
109 if len(dataset.columns) < 2: 109 if len(dataset.columns) < 2:
110 sys.exit('Execution aborted: wrong format of ' + name + '\n') 110 sys.exit('Execution aborted: wrong format of ' + name + '\n')
111 return dataset 111 return dataset
112 112
113 113
114 def apply_ras_bounds(model, ras_row, rxns_ids, mediumRxns_ids): 114 def apply_ras_bounds(model, ras_row):
115 """ 115 """
116 Adjust the bounds of reactions in the model based on RAS values. 116 Adjust the bounds of reactions in the model based on RAS values.
117 117
118 Args: 118 Args:
119 model (cobra.Model): The metabolic model to be modified. 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. 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: 121 Returns:
124 None 122 None
125 """ 123 """
126 for reaction in ras_row.index: 124 for reaction in ras_row.index:
127 scaling_factor = ras_row[reaction] 125 scaling_factor = ras_row[reaction]
136 if upper_bound!=0 and lower_bound!=0: 134 if upper_bound!=0 and lower_bound!=0:
137 model.reactions.get_by_id(reaction).lower_bound=valMin 135 model.reactions.get_by_id(reaction).lower_bound=valMin
138 model.reactions.get_by_id(reaction).upper_bound=valMax 136 model.reactions.get_by_id(reaction).upper_bound=valMax
139 pass 137 pass
140 138
141 def process_ras_cell(cellName, ras_row, model, rxns_ids, mediumRxns_ids, output_folder): 139 def process_ras_cell(cellName, ras_row, model, rxns_ids, output_folder):
142 """ 140 """
143 Process a single RAS cell, apply bounds, and save the bounds to a CSV file. 141 Process a single RAS cell, apply bounds, and save the bounds to a CSV file.
144 142
145 Args: 143 Args:
146 cellName (str): The name of the RAS cell (used for naming the output file). 144 cellName (str): The name of the RAS cell (used for naming the output file).
147 ras_row (pd.Series): A row from a RAS DataFrame containing scaling factors for reaction bounds. 145 ras_row (pd.Series): A row from a RAS DataFrame containing scaling factors for reaction bounds.
148 model (cobra.Model): The metabolic model to be modified. 146 model (cobra.Model): The metabolic model to be modified.
149 rxns_ids (list of str): List of reaction IDs to which the scaling factors will be applied. 147 rxns_ids (list of str): List of reaction IDs to which the scaling factors will be applied.
150 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.
151 output_folder (str): Folder path where the output CSV file will be saved. 148 output_folder (str): Folder path where the output CSV file will be saved.
152 149
153 Returns: 150 Returns:
154 None 151 None
155 """ 152 """
156 model_new = model.copy() 153 model_new = model.copy()
157 apply_ras_bounds(model_new, ras_row, rxns_ids, mediumRxns_ids) 154 apply_ras_bounds(model_new, ras_row)
158 bounds = pd.DataFrame([(rxn.lower_bound, rxn.upper_bound) for rxn in model_new.reactions], index=rxns_ids, columns=["lower_bound", "upper_bound"]) 155 bounds = pd.DataFrame([(rxn.lower_bound, rxn.upper_bound) for rxn in model_new.reactions], index=rxns_ids, columns=["lower_bound", "upper_bound"])
159 bounds.to_csv(output_folder + cellName + ".csv", sep='\t', index=True) 156 bounds.to_csv(output_folder + cellName + ".csv", sep='\t', index=True)
160 pass 157 pass
161 158
162 def generate_bounds(model: cobra.Model, medium: dict, ras=None, output_folder='output/') -> pd.DataFrame: 159 def generate_bounds(model: cobra.Model, medium: dict, ras=None, output_folder='output/') -> pd.DataFrame:
172 Returns: 169 Returns:
173 pd.DataFrame: DataFrame containing the bounds of reactions in the model. 170 pd.DataFrame: DataFrame containing the bounds of reactions in the model.
174 """ 171 """
175 rxns_ids = [rxn.id for rxn in model.reactions] 172 rxns_ids = [rxn.id for rxn in model.reactions]
176 173
177 #metto a zero medium engro2 174 # Set all reactions to zero in the medium
178 for rxn_id, _ in model.medium.items(): 175 for rxn_id, _ in model.medium.items():
179 model.reactions.get_by_id(rxn_id).lower_bound = float(0.0) 176 model.reactions.get_by_id(rxn_id).lower_bound = float(0.0)
180 177
181 # Set medium conditions 178 # Set medium conditions
182 for reaction, value in medium.items(): 179 for reaction, value in medium.items():
183 if value is not None: 180 if value is not None:
184 ## SOLO ENGRO2 181 model.reactions.get_by_id(reaction).lower_bound = -float(value)
185 if(reaction != "EX_thbpt_e" and reaction != "EX_lac__L_e"): 182
186 model.reactions.get_by_id(reaction).lower_bound = -float(value)
187 if(reaction == "EX_lac__L_e"):
188 model.reactions.get_by_id(reaction).lower_bound = float(0.0)
189 183
190 mediumRxns_ids = medium.keys() 184 mediumRxns_ids = medium.keys()
191 185
192 # Perform Flux Variability Analysis (FVA) on this medium 186 # Perform Flux Variability Analysis (FVA) on this medium
193 df_FVA = cobra.flux_analysis.flux_variability_analysis(model, fraction_of_optimum=0, processes=1).round(8) 187 df_FVA = cobra.flux_analysis.flux_variability_analysis(model, fraction_of_optimum=0, processes=1).round(8)
198 model.reactions.get_by_id(reaction).upper_bound = float(df_FVA.loc[reaction, "maximum"]) 192 model.reactions.get_by_id(reaction).upper_bound = float(df_FVA.loc[reaction, "maximum"])
199 193
200 if ras is not None: 194 if ras is not None:
201 #Parallel(n_jobs=cpu_count())(delayed(process_ras_cell)(cellName, ras_row, model, rxns_ids, output_folder) for cellName, ras_row in ras.iterrows()) 195 #Parallel(n_jobs=cpu_count())(delayed(process_ras_cell)(cellName, ras_row, model, rxns_ids, output_folder) for cellName, ras_row in ras.iterrows())
202 for cellName, ras_row in ras.iterrows(): 196 for cellName, ras_row in ras.iterrows():
203 process_ras_cell(cellName, ras_row, model, rxns_ids, mediumRxns_ids, output_folder) 197 process_ras_cell(cellName, ras_row, model, rxns_ids, output_folder)
204 break #just one cell for testing 198 break #just one cell for testing
205 else: 199 else:
206 model_new = model.copy() 200 model_new = model.copy()
207 apply_ras_bounds(model_new, pd.Series([1]*len(rxns_ids), index=rxns_ids), rxns_ids, mediumRxns_ids) 201 apply_ras_bounds(model_new, pd.Series([1]*len(rxns_ids), index=rxns_ids))
208 bounds = pd.DataFrame([(rxn.lower_bound, rxn.upper_bound) for rxn in model_new.reactions], index=rxns_ids, columns=["lower_bound", "upper_bound"]) 202 bounds = pd.DataFrame([(rxn.lower_bound, rxn.upper_bound) for rxn in model_new.reactions], index=rxns_ids, columns=["lower_bound", "upper_bound"])
209 bounds.to_csv(output_folder + "bounds.csv", sep='\t', index=True) 203 bounds.to_csv(output_folder + "bounds.csv", sep='\t', index=True)
210 pass 204 pass
211 205
212 206
250 # Concatenate all ras DataFrames into a single DataFrame 244 # Concatenate all ras DataFrames into a single DataFrame
251 ras_combined = pd.concat(ras_list, axis=0) 245 ras_combined = pd.concat(ras_list, axis=0)
252 # Normalize the RAS values by max RAS 246 # Normalize the RAS values by max RAS
253 ras_combined = ras_combined.div(ras_combined.max(axis=0)) 247 ras_combined = ras_combined.div(ras_combined.max(axis=0))
254 ras_combined.dropna(axis=1, how='all', inplace=True) 248 ras_combined.dropna(axis=1, how='all', inplace=True)
255 #ras_combined = ras_combined.fillna(0)
256 #il ras c'è per tutti o non c'è per nessuno
257 249
258 250
259 251
260 model_type :utils.Model = ARGS.model_selector 252 model_type :utils.Model = ARGS.model_selector
261 if model_type is utils.Model.Custom: 253 if model_type is utils.Model.Custom: