#!/usr/bin/python from cbraPytool import * # Define minimal media by setting maximum uptake rate for exchange metabolites (nutrients) # Minimal media uptakeBound_mm = {'M_o2_b': 6.3, 'M_nh4_b': 100, # Ammonium (nitrogen source) 'M_so4_b': 100, # sulfate 'M_pi_b': 0.89, # phospate 'M_glc_D_b': 22.6, # D-glucose 'M_his_L_b': 0.082, # histidin 'M_leu_L_b': 0.4, # Leu 'M_met_L_b': 0.0468, # Met 'M_k_b': 4.44, # Potassium "M_na1_b": 0.75, # Sodium "M_btn_b": 0.00000142, #Biotin "M_chol_b": 0.000092, "M_inost_b": 0.00193, "M_pnto_R_b": 0.0002, "M_ribflv_b": 0.00063, "M_ura_b": 0.4 } # Rich (YPD) media uptakeBound_ypd = { 'M_o2_b': 6.3, 'M_nh4_b': 100, 'M_so4_b': 100, 'M_pi_b': 0.89, 'M_glc_D_b': 22.6, 'M_ala_L_b': 0.519, "M_arg_L_b": 0.193, "M_asp_L_b": 0.274, "M_cys_L_b": 0.0192, "M_glu_L_b": 0.479, "M_gly_b": 0.934, 'M_his_L_b': 0.31, "M_ile_L_b": 0.0952, 'M_leu_L_b': 1.66, "M_lys_L_b": 0.167, "M_met_L_b": 0.0468, "M_phe_L_b": 0.0758, "M_pro_L_b": 0.357, "M_ser_L_b": 0.166, "M_thr_L_b": 0.112, "M_trp_L_b": 0.0207, "M_tyr_L_b": 0.0279, "M_val_L_b": 0.148, 'M_k_b': 1.87, "M_na1_b": 4.43, "M_btn_b": 0.0000275, "M_chol_b": 0.00437, "M_ribflv_b": 0.00063, "M_thmpp_b": 0.0032, "M_inost_b": 0.07, "M_thymd_b": 0.00709, "M_pnto_R_b": 0.0002, "M_ura_b": 4} print "-- Define wild-type auxotrophic marker deletion" KO0 = 'YEL021W:YLR303W:YCL018W:YOR202W' print "-- Get a list of genes of interest" lstKO = ['', 'YDL052C',#2.3.1.51,1-acylglycerol-3-phosphate O-acyltransferase, 'R_AGAT_SC': existing annotation!!! 'YDL168W',# 1.1.1.103,L-threonine 3-dehydrogenase, No reaction, only L_threonine_deaminase, "R_THRD_L" has 2 associated genes already 'YER152C',#2.6.1.39, 2-aminoadipate transaminase,'R_AATA', not exactly the same reaction 'YER183C',#6.3.3.2, 5,10-methenyltetrahydrofolate synthetase, 'R_FTHFCL', existing annotation!!! ; R_FTHFCLm, the gene annotation for mitochondrion is missing, does it stay in mitochondiron???? 'YGL202W',#2.6.1.39, 2-aminoadipate transaminase,'R_AATA', not exactly the same reaction 'YGR248W',#3.5.99.6,glucosamine-6-phosphate deaminase, 'R_G6PDA' 'YHR163W',#3.5.99.6,glucosamine-6-phosphate deaminase, 'R_G6PDA' 'YIL033C',#3.5.1.2,glutaminase, R_GLUN 'YJL045W',#1.4.3.16,L-aspartate oxidase, No such reaction 'YJL060W',#2.6.1.39, 2-aminoadipate transaminase,'R_AATA' 'YJL218W',#2.3.1.30, serine O-acetyltransferase,'R_SERATi' 'YLL060C',#5.2.1.2, maleylacetoacetate isomerase, 'R_MACACI' 'YLR017W',# 2.4.2.1, 'R_PUNP1'??? no 5-methylthio though, adding isoenzyme to R_PUNP1?? 'YLR070C',# 1.1.1.103,L-threonine 3-dehydrogenase, No reaction, only L_threonine_deaminase, "R_THRD_L" has 2 associated genes already 'YLR209C',#2.4.2.1, purine-nucleoside phosphorylase, existing annotation but in different forms of reaction: with added metabolites e.g. R_PUNP2: "R_purine_nucleoside_phosphorylase__Deoxyadenosine_", R_purine_nucleoside_phosphorylase__Xanthosine_, R_PUNP1: R_purine_nucleoside_phosphorylase__Adenosine_, R_PUNP1m: R_purine_nucleoside_phosphorylase__Adenosine___mitochondrial, R_PUNP3:R_purine_nucleoside_phosphorylase__Guanosine_, R_PUNP3m:R_purine_nucleoside_phosphorylase__Guanosine___mitochondrial, R_PUNP4:R_purine_nucleoside_phosphorylase__Deoxyguanosine_, R_PUNP5:R_purine_nucleoside_phosphorylase__Inosine_, R_PUNP6:R_purine_nucleoside_phosphorylase__Deoxyinosine_, R_PUNP7:R_purine_nucleoside_phosphorylase__Xanthosine_, 'YMR020W', # 1.5.3.11, poylamine_oxidase, 'R_POLYAO2' 'YNR027W', # 2.7.1.35, pyridoxal kinase,R_PYDXK 'YNR034W', # 3.5.99.6, glucosamine-6-phosphate deaminase, 'R_G6PDA' 'YNR073C', # 1.1.1.17, mannitol-1-phosphate 5-dehydrogenase, "R_M1PD" 'YPR121W'] # 2.7.1.35, pyridoxal kinase,R_PYDXK #lstKO0 = lstKO[1:] ######## print "-- Define the putative annotation from Adam's experiment" # For multiple gene anototated reaction, it is assumed that all enzymes ecoded by the genes # should be present in order to catalyze the reaction. # All enzymes involved in multiple genes are treated as enzyme complexes here as can be seen # genes are separated by ':' # To use iso-enzyme treatment for the annotated genes, use '|' to separate them instead of ':' geneAnnot = { 'R_AGAT_SC':'YDL052C', 'R_M1PD':'YNR073C', 'R_POLYAO2':'YMR020W', 'R_SERATi':'YJL218W', 'R_AATA':'YER152C:YGL202W:YJL060W', 'R_GLUN':'YIL033C', 'R_G6PDA':'YHR163W:YNR034W:YGR248', 'R_MACACI':'YLL060C', 'R_PYDXK':'YNR027W:YPR121W', 'R_FTHFCLm':'YER183C', 'R_PUNP1':'YLR209C:YLR017W', # or original 'YLR209C' } def main(): # Reading model from xml file print "-- Create a new empty constraint-based metabolic model" cbm0 = CbModel() print "-- Read SBML model into CbModel format in python" cbm0.readSbmlCbModel(sbmlFile="Sc_iND750_GlcMM.xml") print "-- Make a copy of the original model cbm0" cbm1 = cbm0.copy() print "-- wild-type simulation (without any auxotrophic marker deletion) under minimal media" print "-- set up the maximum uptake rate for nutrients from the media" cbm1.setUptakeBound(uptakeBound_mm) print "-- running standard FBA simulation for wild-type and store solution in sol_wt" sol_wt = cbm1.optCbModel(solver='lpsolve_orig') print "-- wildtype with auxotrophic marker deletion" cbm1 = cbm1.deleteGenes(genes = KO0) print "-- running standard FBA simulation for wild-type and store" sol_wtax=cbm1.optCbModel(solver='lpsolve_orig') nR = len(cbm1.lstR) # Modify the relevant gene annotations for reuse study # the new annotation list (lstEnz) will be used for deletion simulation lstEnz = deepcopy(cbm0.lstEnz) lstEnz0 = cbm0.lstEnz for r1 in geneAnnot: idx=arange(nR)[array(cbm0.lstR)==r1] lstEnz[idx] = geneAnnot[r1] print 'Annotation for reaction %s: %s -> %s' %(r1, cbm1.dict['O'][r1]['Enz'], geneAnnot[r1]) cbm1.dict['O'][r1]['Enz'] = geneAnnot[r1] cbm1.lstEnz = lstEnz # Flux range analysis under different media and different constraints # FVA for reactions involved. lstRD = geneAnnot.keys() nR = len(cbm0.lstR) lstv = [arange(nR).compress(cbm0.lstR==k)[0] for k in lstRD] R_growth = cbm1.lstR[cbm1.objC == 1][0] # get ID for biomass formation reaction for media in ['sd', 'ypd']: print media # Set up uptake contraints if media[0:2] == 'sd': uptakeBound = uptakeBound_mm elif media[0:3] == 'ypd': # YPD medium uptakeBound = uptakeBound_ypd cbm1.setUptakeBound(uptakeBound = uptakeBound) print"-- Structural constraints " [lstvMin, lstvMax, lstObjMin, lstObjMax, lst_sol_min, lst_sol_max, lstv] = cbm1.fva(lstv=lstv, fileOut='fva_reuse_'+media+'StructuralConstr.csv', fixedObjC=False, constrR=None) print"-- Structural contraints and minimal grow th rate of 0.05" [lstvMin1, lstvMax1, lstObjMin1, lstObjMax1, lst_sol_min, lst_sol_max, lstv] = cbm1.fva(lstv=lstv, fileOut='fva_'+media+'MinimalGrowth.csv', fixedObjC=False, constrR={R_growth:{'lb':0.05}}) ###################### # Prediction for growth rate change for deletion mutant compared to wild-type. ###################### method = 'fba' #method = 'moma' array_biomass = [] array_sol = [] for media in ['sd', 'ypd']: print 'Media: %s' %(media) # Set up uptake contraints if media[0:2] == 'sd': uptakeBound = uptakeBound_mm elif media[0:3] == 'ypd': uptakeBound = uptakeBound_ypd cbm1.setUptakeBound(uptakeBound = uptakeBound) print len(uptakeBound), '!!!' print 'Start testing ...' [lstGr1, lstSol1, null] = cbm1.simGeneDeletion(lstGene=lstKO, fileOut='deletion_'+method+'_ind750_'+media+'.csv', method=method) array_biomass.append(lstGr1) array_sol.append(lstSol1) if __name__ == '__main__': main()