|
93
|
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
|
|
98
|
57 parser.add_argument('-rn', '--name',
|
|
94
|
58 type=str,
|
|
|
59 help = 'ras class names')
|
|
|
60
|
|
93
|
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
|
|
94
|
70
|
|
93
|
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
|
|
113
|
114 def apply_ras_bounds(model, ras_row, rxns_ids, mediumRxns_ids):
|
|
93
|
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.
|
|
113
|
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.
|
|
93
|
123 Returns:
|
|
|
124 None
|
|
|
125 """
|
|
115
|
126 warning(str(model.reactions.get_by_id("EX_Lcystin_e").lower_bound))
|
|
116
|
127 warning(str(model.reactions.get_by_id("EX_Lcystin_e").upper_bound))
|
|
93
|
128 for reaction in rxns_ids:
|
|
115
|
129 if reaction in ras_row.index: #and reaction not in mediumRxns_ids
|
|
93
|
130 scaling_factor = ras_row[reaction]
|
|
|
131 lower_bound=model.reactions.get_by_id(reaction).lower_bound
|
|
|
132 upper_bound=model.reactions.get_by_id(reaction).upper_bound
|
|
|
133 valMax=float((upper_bound)*scaling_factor)
|
|
|
134 valMin=float((lower_bound)*scaling_factor)
|
|
|
135 if upper_bound!=0 and lower_bound==0:
|
|
|
136 model.reactions.get_by_id(reaction).upper_bound=valMax
|
|
|
137 if upper_bound==0 and lower_bound!=0:
|
|
|
138 model.reactions.get_by_id(reaction).lower_bound=valMin
|
|
|
139 if upper_bound!=0 and lower_bound!=0:
|
|
|
140 model.reactions.get_by_id(reaction).lower_bound=valMin
|
|
|
141 model.reactions.get_by_id(reaction).upper_bound=valMax
|
|
115
|
142 warning(str(model.reactions.get_by_id("EX_Lcystin_e").lower_bound))
|
|
116
|
143 warning(str(model.reactions.get_by_id("EX_Lcystin_e").upper_bound))
|
|
93
|
144 pass
|
|
|
145
|
|
113
|
146 def process_ras_cell(cellName, ras_row, model, rxns_ids, mediumRxns_ids, output_folder):
|
|
93
|
147 """
|
|
|
148 Process a single RAS cell, apply bounds, and save the bounds to a CSV file.
|
|
|
149
|
|
|
150 Args:
|
|
|
151 cellName (str): The name of the RAS cell (used for naming the output file).
|
|
|
152 ras_row (pd.Series): A row from a RAS DataFrame containing scaling factors for reaction bounds.
|
|
|
153 model (cobra.Model): The metabolic model to be modified.
|
|
|
154 rxns_ids (list of str): List of reaction IDs to which the scaling factors will be applied.
|
|
113
|
155 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.
|
|
93
|
156 output_folder (str): Folder path where the output CSV file will be saved.
|
|
|
157
|
|
|
158 Returns:
|
|
|
159 None
|
|
|
160 """
|
|
|
161 model_new = model.copy()
|
|
113
|
162 apply_ras_bounds(model_new, ras_row, rxns_ids, mediumRxns_ids)
|
|
93
|
163 bounds = pd.DataFrame([(rxn.lower_bound, rxn.upper_bound) for rxn in model_new.reactions], index=rxns_ids, columns=["lower_bound", "upper_bound"])
|
|
|
164 bounds.to_csv(output_folder + cellName + ".csv", sep='\t', index=True)
|
|
|
165 pass
|
|
|
166
|
|
|
167 def generate_bounds(model: cobra.Model, medium: dict, ras=None, output_folder='output/') -> pd.DataFrame:
|
|
|
168 """
|
|
|
169 Generate reaction bounds for a metabolic model based on medium conditions and optional RAS adjustments.
|
|
|
170
|
|
|
171 Args:
|
|
|
172 model (cobra.Model): The metabolic model for which bounds will be generated.
|
|
|
173 medium (dict): A dictionary where keys are reaction IDs and values are the medium conditions.
|
|
|
174 ras (pd.DataFrame, optional): RAS pandas dataframe. Defaults to None.
|
|
|
175 output_folder (str, optional): Folder path where output CSV files will be saved. Defaults to 'output/'.
|
|
|
176
|
|
|
177 Returns:
|
|
|
178 pd.DataFrame: DataFrame containing the bounds of reactions in the model.
|
|
|
179 """
|
|
|
180 rxns_ids = [rxn.id for rxn in model.reactions]
|
|
|
181
|
|
|
182 # Set medium conditions
|
|
|
183 for reaction, value in medium.items():
|
|
|
184 if value is not None:
|
|
107
|
185 if(reaction != "EX_thbpt_e" and reaction != "EX_lac__L_e"):
|
|
|
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)
|
|
113
|
189
|
|
|
190 mediumRxns_ids = medium.keys()
|
|
107
|
191
|
|
93
|
192 # Perform Flux Variability Analysis (FVA)
|
|
|
193 df_FVA = cobra.flux_analysis.flux_variability_analysis(model, fraction_of_optimum=0, processes=1).round(8)
|
|
|
194
|
|
|
195 # Set FVA bounds
|
|
|
196 for reaction in rxns_ids:
|
|
102
|
197 model.reactions.get_by_id(reaction).lower_bound = float(df_FVA.loc[reaction, "minimum"])
|
|
|
198 model.reactions.get_by_id(reaction).upper_bound = float(df_FVA.loc[reaction, "maximum"])
|
|
93
|
199
|
|
118
|
200 warning(str(model.reactions.get_by_id("EX_Lcystin_e").lower_bound))
|
|
|
201 warning(str(model.reactions.get_by_id("EX_Lcystin_e").upper_bound))
|
|
|
202
|
|
93
|
203 if ras is not None:
|
|
106
|
204 #Parallel(n_jobs=cpu_count())(delayed(process_ras_cell)(cellName, ras_row, model, rxns_ids, output_folder) for cellName, ras_row in ras.iterrows())
|
|
|
205 for cellName, ras_row in ras.iterrows():
|
|
113
|
206 process_ras_cell(cellName, ras_row, model, rxns_ids, mediumRxns_ids, output_folder)
|
|
107
|
207 break #just one cell for testing
|
|
93
|
208 else:
|
|
|
209 model_new = model.copy()
|
|
113
|
210 apply_ras_bounds(model_new, pd.Series([1]*len(rxns_ids), index=rxns_ids), rxns_ids, mediumRxns_ids)
|
|
93
|
211 bounds = pd.DataFrame([(rxn.lower_bound, rxn.upper_bound) for rxn in model_new.reactions], index=rxns_ids, columns=["lower_bound", "upper_bound"])
|
|
|
212 bounds.to_csv(output_folder + "bounds.csv", sep='\t', index=True)
|
|
|
213 pass
|
|
|
214
|
|
|
215
|
|
|
216
|
|
|
217 ############################# main ###########################################
|
|
|
218 def main() -> None:
|
|
|
219 """
|
|
|
220 Initializes everything and sets the program in motion based on the fronted input arguments.
|
|
|
221
|
|
|
222 Returns:
|
|
|
223 None
|
|
|
224 """
|
|
|
225 if not os.path.exists('ras_to_bounds'):
|
|
|
226 os.makedirs('ras_to_bounds')
|
|
|
227
|
|
|
228
|
|
|
229 global ARGS
|
|
|
230 ARGS = process_args(sys.argv)
|
|
|
231
|
|
|
232 ARGS.output_folder = 'ras_to_bounds/'
|
|
|
233
|
|
|
234 if(ARGS.ras_selector == True):
|
|
|
235 ras_file_list = ARGS.input_ras.split(",")
|
|
98
|
236 ras_file_names = ARGS.name.split(",")
|
|
94
|
237 ras_class_names = []
|
|
|
238 for file in ras_file_names:
|
|
|
239 ras_class_names.append(file.split(".")[0])
|
|
93
|
240 ras_list = []
|
|
|
241 class_assignments = pd.DataFrame(columns=["Patient_ID", "Class"])
|
|
|
242 for ras_matrix, ras_class_name in zip(ras_file_list, ras_class_names):
|
|
|
243 ras = read_dataset(ras_matrix, "ras dataset")
|
|
|
244 ras.replace("None", None, inplace=True)
|
|
|
245 ras.set_index("Reactions", drop=True, inplace=True)
|
|
|
246 ras = ras.T
|
|
|
247 ras = ras.astype(float)
|
|
|
248 ras_list.append(ras)
|
|
|
249 for patient_id in ras.index:
|
|
94
|
250 class_assignments.loc[class_assignments.shape[0]] = [patient_id, ras_class_name]
|
|
|
251
|
|
93
|
252
|
|
|
253 # Concatenate all ras DataFrames into a single DataFrame
|
|
94
|
254 ras_combined = pd.concat(ras_list, axis=0)
|
|
93
|
255 # Normalize the RAS values by max RAS
|
|
|
256 ras_combined = ras_combined.div(ras_combined.max(axis=0))
|
|
|
257 ras_combined = ras_combined.fillna(0)
|
|
|
258
|
|
|
259
|
|
|
260
|
|
|
261 model_type :utils.Model = ARGS.model_selector
|
|
|
262 if model_type is utils.Model.Custom:
|
|
|
263 model = model_type.getCOBRAmodel(customPath = utils.FilePath.fromStrPath(ARGS.model), customExtension = utils.FilePath.fromStrPath(ARGS.model_name).ext)
|
|
|
264 else:
|
|
|
265 model = model_type.getCOBRAmodel(toolDir=ARGS.tool_dir)
|
|
|
266
|
|
|
267 if(ARGS.medium_selector == "Custom"):
|
|
|
268 medium = read_dataset(ARGS.medium, "medium dataset")
|
|
|
269 medium.set_index(medium.columns[0], inplace=True)
|
|
|
270 medium = medium.astype(float)
|
|
|
271 medium = medium[medium.columns[0]].to_dict()
|
|
|
272 else:
|
|
|
273 df_mediums = pd.read_csv(ARGS.tool_dir + "/local/medium/medium.csv", index_col = 0)
|
|
|
274 ARGS.medium_selector = ARGS.medium_selector.replace("_", " ")
|
|
|
275 medium = df_mediums[[ARGS.medium_selector]]
|
|
|
276 medium = medium[ARGS.medium_selector].to_dict()
|
|
|
277
|
|
|
278 if(ARGS.ras_selector == True):
|
|
|
279 generate_bounds(model, medium, ras = ras_combined, output_folder=ARGS.output_folder)
|
|
94
|
280 class_assignments.to_csv(ARGS.cell_class, sep = '\t', index = False)
|
|
93
|
281 else:
|
|
|
282 generate_bounds(model, medium, output_folder=ARGS.output_folder)
|
|
|
283
|
|
|
284 pass
|
|
|
285
|
|
|
286 ##############################################################################
|
|
|
287 if __name__ == "__main__":
|
|
|
288 main() |
