<jupyter_start><jupyter_text>Figure 1 from Zhao 2019![figure1](figure1-Zhao2019.png)<jupyter_code>import numpy as np
np.random.seed(0)

import pandas as pd

def convert_timepoints_to_flux( timepoints_ref ):
    timepoints = timepoints_ref.copy()
    return pd.DataFrame( timepoints.loc[:, timepoints.columns[1:]].values - timepoints.loc[:, timepoints.columns[:-1]].values,
                         index = timepoints.index, columns = timepoints.columns[1:])

%matplotlib inline<jupyter_output><empty_output><jupyter_text>Here we construct the model using the `cobrapy` package, to get the stochiometric matrix from the list of reactions  Read Table 3 from Zhao 2019<jupyter_code>table3 = pd.read_csv('Table3-Zhao2019.tab', sep='\t', index_col=['metabolite','mean-or-std'])
table3<jupyter_output><empty_output><jupyter_text>Average daily metabolite flux across all strains<jupyter_code>idx = pd.IndexSlice                                                         # This allows us to select into a multi-index row index.
table3_mean_conc = table3.loc[ idx[ :, 'mean' ], : ]                             # return a table that only includes the mean of each metabolite, not the std
table3_mean_conc.index = table3_mean_conc.index.droplevel('mean-or-std')              # remove the mean-or-std row index.
table3_mean_conc = table3_mean_conc.drop('Unnamed: 0',axis=1)   # drop the first column
table3_mean_conc.index.name = 'Metabolite $g/L$'
table3_flux = table3_mean_conc.copy()
table3_flux = table3_flux.fillna(0)/4.0                                     # replace NaN with 0 and divide everything by 4 (since we want flux/day and concentrations were taken after 96 hours)
table3_flux.index.name = 'Metabolite $g/L/day$'                             # Provide informative units for the row index name
table3_flux.to_csv('table3_zhao2019.strains.grams_per_liter_per_day.tab', sep='\t')  # Print out to file
table3_flux
table3_mean_conc<jupyter_output><empty_output><jupyter_text>Inferred enzyme overexpression levels from Table 3<jupyter_code>enzyme_table3 = pd.read_csv('enzyme-overexpression-Zhao2019.tab', sep='\t',index_col=0)
enzyme_table3.index = ['CAD'] + list(enzyme_table3.index)[1:]
enzyme_table3.loc['CAD'] = 1
enzyme_table3.loc['CAD', 'Po1f'] = 0
enzyme_table3.loc['IDH'] = 0.5
enzyme_table3.loc['MTT','Po1f'] = 0
enzyme_table3.loc['MFS'] = 0.5
enzyme_table3.loc['MFS', ['PC-MFS','PCM-MFS', 'PCM-ACO-MFS']]= 1
enzyme_table3.to_csv('enzyme-overexpression-Zhao2019.tab', sep='\t')
enzyme_table3<jupyter_output><empty_output><jupyter_text>Table 4: Timeseries concentration of metabolites from PCM-ACO-MFS<jupyter_code>table4 =pd.read_csv('Table4-Zhao2019.tab', sep='\t', index_col=['metabolite','mean-or-std'])
table4
%matplotlib inline
import plotnine as p9
from plotnine import ggplot, aes, geom_point, geom_line
idx = pd.IndexSlice
itacon_timeseries = table4.loc[idx['Itaconic acid','mean'],:]
ts = itacon_timeseries.index.insert(2, '96 h')
itacon_timeseries['96 h'] = 5.21
itacon_timeseries = itacon_timeseries.reindex(ts)
itacon_timeseries.index = [int(i[:-1].strip()) for i in itacon_timeseries.index]
itacon_timeseries.plot()
itacon_timeseries.to_frame()<jupyter_output><empty_output><jupyter_text>Timeseries flux of metabolites for PCM-ACO-MFSTo convert metabolite concentrations in $g/L$ to metabolite fluxes in $g/L/day$ we must weight the difference in concentration by the number of hours between measurements$$v_{t_n} = \left(c_{t_n} -c_{t_{n-1}}\right)\left(\frac{24}{t_{n} - t_{n-1}}\right)$$where $n$ ranges from timepoint 0 to timepoint $N$, $t_n$ is the number of hours at timepoint $n$, $c_{t_n}$ is the concentration at timepoint $n$, and $v_{t_n}$ is the flux at timepoint $n$.  If the number of concentration timepoints range from $0..N$, the number of flux timepoints range from $1..N$.<jupyter_code>table4_mean_conc = table4.loc[ idx[ :, 'mean' ], : ]                      # Generate a table that only contains mean concentrations.
table4_mean_conc.index = table4_mean_conc.index.droplevel('mean-or-std')  # drop the mean-or-std row index
table4_mean_conc.index.name = 'Metabolite $g/L$'
timepoints = np.array( [ int( col[:-1].strip() )                          # Generate timepoints from column index
                        for col in table4_mean_conc.columns] ) 
table4_flux = convert_timepoints_to_flux( table4_mean_conc.copy() )\
                    .multiply( 24/( timepoints[1:] - timepoints[:-1] ))   # apply the formula converting concentrations and timepoints to fluxes
table4_flux.index.name= 'Metabolite $g/L/day$'                       # Provide informative units for the row index name
table4_flux = table4_flux.fillna(0)                             # convert NaN to 0
table4_flux.to_csv('table4_zhao2019.PCM-ACO-MFS.timeseries.grams_per_liter_per_day.tab', 
                         sep='\t')                                        # write out the dataframe to a file.
table4_flux
table4_mean_conc<jupyter_output><empty_output><jupyter_text>Inferred enzyme overexpression from table 4<jupyter_code>enzyme_table4 = pd.DataFrame({time_point: enzyme_table3['PCM-ACO-MFS'] 
                              for time_point in table4.columns})
enzyme_table4<jupyter_output><empty_output><jupyter_text>Figure 3![figure3](figure3-Zhao2019.png)  Itaconate acid production in *Y. lipolytica* strains expressing the *cadA* gene.<jupyter_code>figure3 = pd.read_excel('Flux_data_Zhao2019.xlsx', 
                         sheet_name = 'Figure3',index_col=0)
figure3.index.name = 'Metabolite $g/L$'
figure3.to_csv('figure3.cadA.timeseries.grams_per_liter.tab', sep='\t')
figure3<jupyter_output><empty_output><jupyter_text>Daily flux of *cadA* containing strains computed from average difference of metabolite concentrations in g/L/day<jupyter_code>figure3_flux = convert_timepoints_to_flux( figure3 )
figure3_flux.index.name = 'Metabolite $g/L/day$'
figure3_flux.to_csv('figure3.cadA.timeseries.grams_per_liter_per_day.tab', sep='\t')
figure3_flux<jupyter_output><empty_output><jupyter_text>Figure 4![figure4](figure4-Zhao2019.png)<jupyter_code>figure4 = {}
for sheet_name in 'abcd':
    figure4[sheet_name] = pd.read_excel('Flux_data_Zhao2019.xlsx', 
                                         sheet_name = 'Figure4{}'.format( sheet_name ),
                                       index_col=0)
    figure4[sheet_name].index.name = 'Metabolite $g/L$'<jupyter_output><empty_output><jupyter_text>Figure 4a: Comparison of different genetic modifications in *Y. lipolytica* PC-7 across 4 days<jupyter_code>figure4['a'].to_csv('figure4a_zhao2019.PC.children.grams_per_liter.tab', sep='\t')
figure4['a']<jupyter_output><empty_output><jupyter_text>Flux per day averaged across 4 days for  metabolites of PC strains in g/L/day<jupyter_code>figure4a_flux = figure4['a']/4
figure4a_flux.index.name = 'Metabolite $g/L/day$'
figure4a_flux.to_csv('figure4a_zhao2019.PC.children.grams_per_liter_per_day.tab', sep='\t')
figure4a_flux<jupyter_output><empty_output><jupyter_text>Figure 4b: Shake flask fermentation results of the PC-MTT<jupyter_code>figure4['b'].to_csv('figure4b_zhao2019.PC-MTT.timeseries.grams_per_liter.tab', sep='\t')
figure4['b']<jupyter_output><empty_output><jupyter_text>Daily flux of PC-MTT strain computed from average difference of metabolite concentrations in $g/L/day$<jupyter_code>figure4b_flux = convert_timepoints_to_flux( figure4['b'])
figure4b_flux.index.name = 'Metabolite $g/L/day$'
figure4b_flux.to_csv('figure4b_zhao2019.PC-MTT.timeseries.grams_per_liter_per_day.tab', sep='\t')
figure4b_flux<jupyter_output><empty_output><jupyter_text>Figure 4c: Comparison of different genetic modifications in *Y. lipolytica* PC-MTT<jupyter_code>figure4['c'].to_csv('figure4c_zhao2019.PC-MTT.children.grams_per_liter.tab',sep='\t')
figure4['c']<jupyter_output><empty_output><jupyter_text>Flux per day averaged across 4 days for metabolites of children of PC-MTT parent strain in $g/L/day$<jupyter_code>figure4c_flux = figure4['c']/4
figure4c_flux.index.name = 'Metabolite $g/L/day$'
figure4c_flux.to_csv('figure4c_zhao2019.PC-MTT.children.grams_per_liter_per_day.tab', sep='\t')
figure4c_flux<jupyter_output><empty_output><jupyter_text>Figure 4d: Shake flask fermentation results of PCM-ACO-MFS<jupyter_code>figure4['d'].to_csv('figure4d-zhao2019.PCM-ACO-MFS.timeseries.grams_per_liter.tab', sep='\t')
figure4['d']<jupyter_output><empty_output><jupyter_text>Daily flux of PCM-ACO-MFS strain computed from average difference of metabolite concentrations in $g/L/day$<jupyter_code>figure4d_flux = convert_timepoints_to_flux( figure4['d'] )
figure4d_flux.index.name = 'Metabolite $g/L/day$'
figure4d_flux.to_csv('figure4d-zhao2019.PCM-ACO-MFS.timeseries.grams_per_liter_per_day.tab', sep='\t')
figure4d_flux<jupyter_output><empty_output><jupyter_text>Map measured metabolite and protein IDs to model<jupyter_code>grams_per_mol = {'Itaconic acid':128.08,
                 'Citric acid': 189.101,
                 'alpha-ketoglutaric acid':144.084,
                 'Glucose': 180.157,
                 'Pyruvate': 87.055,
                 'cis-aconitic acid': 170.086,
                 'trans-aconitic acid': 170.086}<jupyter_output><empty_output><jupyter_text>Map measured metabolites to model metabolites<jupyter_code>measured_metabolites = table4_mean_conc.join(table3_mean_conc).T
measured_metabolites<jupyter_output><empty_output><jupyter_text>According to Table3 of [Zhao 2019], all metabolite measurements were extracelluar.According to [Rafelski 2012](https://f1000.com/fulltext/doi/10.1126/science.1225720), the ratio of cytosol to mitochondrion in yeast is 100:1, so total concentration is mostly cytosolic, and estimates of mitochondrial concentrations would be unreliable.<NAME>., <NAME>., <NAME>., <NAME>., <NAME>., <NAME>., <NAME>., <NAME>., <NAME>. and <NAME>. (2012). Mitochondrial network size scaling in budding yeast. Science 338(6108), 822–824.<jupyter_code>metabolite_map = {
 'pyr_e': 'Pyruvate',
 'acon_C_e': 'cis-aconitic acid',
 'akg_e': 'alpha-ketoglutaric acid',
 'acon_T_e':    'trans-aconitic acid',
 'cit_e': 'Citric acid',
 'itacon_e': 'Itaconic acid'
}
model_metabolites = pd.DataFrame()
for metabolite in metabolite_map:
    model_metabolites[metabolite] = measured_metabolites[metabolite_map[metabolite]]/grams_per_mol[metabolite_map[metabolite]]
model_metabolites.index.name = 'Metabolite $mol/L$'
model_metabolites = model_metabolites.sort_index(axis=1)
model_metabolites.drop('Po1f').drop('20 h').to_csv('glycolysis_itaconate_pwy_metabolites_zhao2019.csv')
model_metabolites<jupyter_output><empty_output><jupyter_text>Map measured enzymes to model reactions<jupyter_code>measured_enzymes = (  set( enzyme_table3.index ) 
                    | set( enzyme_table4.index ) )
measured_enzymes<jupyter_output><empty_output><jupyter_text>Model reaction to Measured Enzyme mapIf multiple enzymes catalyze a reaction, <NAME>. John averaged their expression.We may choose a more sophisticated approach, where we take the min of complex subunits and the sum of isozome expression.$$\begin{eqnarray}  n A  + m B & \rightarrow &  AB \\Expression(complex)  & = &  \min(A/n, B/m)  \\ Expression(isozyme) &=& A + B \\ Expression(multifunctional enzyme) = \sum_{r} \end{eqnarray}$$Since ACO has a mitochondrial location signal, we assume the measured ACO is predominantly in the mitochondrion.<jupyter_code>enzyme_map= {'ACONTa':'ACO',
             'ACONTb':'ACO',
             'MTT':'MTT',
             'ITACON':'CAD',
             'ITACONtr':'MFS',
             'CSm': 'CIT',
             'ICDHym': 'IDH',
             'PPCK': 'PCK',
             'PC': 'PYC'
            }
measured_enzymes = enzyme_table4.join(enzyme_table3).T
#measured_enzymes['CAD'] = 1
#measured_enzymes.loc['Po1f','CAD'] = 0
measured_enzymes
model_enzymes = pd.DataFrame()
for enzyme in enzyme_map:
    model_enzymes[enzyme] = measured_enzymes[enzyme_map[enzyme]]
model_enzymes.index.name = 'Enzyme fold change'
model_enzymes = model_enzymes.sort_index(axis=1)
model_enzymes = model_enzymes.drop('Po1f').drop('20 h')
model_enzymes.to_csv('glycolysis_itaconate_pwy_enzymes_zhao2019.csv')
model_enzymes<jupyter_output><empty_output><jupyter_text>Map measured fluxes to model reactions<jupyter_code>figure4d_flux
measured_flux = table4_flux.join(table3_flux).T
measured_flux['Glucose'] = 21
measured_flux.loc['44 h','Glucose'] = 25
measured_flux.loc[['102 h', '201 h', '308 h', '386 h', '460 h'],'Glucose'] = 0
measured_flux
def model_flux_from_measured_flux( measured_flux, flux_map, grams_per_mol, units =  'Metabolite flux $mol/L/day$' ):
    model_flux = pd.DataFrame()
    for rxn in flux_map:
        model_flux[rxn] = measured_flux[flux_map[rxn]]/grams_per_mol[flux_map[rxn]]
    model_flux = model_flux.sort_index(axis=1)
    model_flux.index.name = units
    return model_flux
flux_map= {
 'EX_itacon_e': 'Itaconic acid',
 'EX_cit_e':    'Citric acid',
 'EX_akg_e':   'alpha-ketoglutaric acid',
 'EX_acon_C_e': 'cis-aconitic acid',
 'EX_glc_D_e':  'Glucose' , 
 'EX_pyr_e':    'Pyruvate',
 'EX_acon_T_e': 'trans-aconitic acid'
}

grams_per_mol = {'Itaconic acid':128.08,
                 'Citric acid': 189.101,
                 'alpha-ketoglutaric acid':144.084,
                 'cis-aconitic acid': 171.09,
                 'Glucose': 180.157,
                 'Pyruvate': 87.055,
                 'trans-aconitic acid': 171.09,
}

model_flux = model_flux_from_measured_flux( measured_flux, flux_map, grams_per_mol )
                    

model_flux.drop('Po1f').to_csv('glycolysis_itaconate_pwy_flux_zhao2019.csv')
model_flux<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.1236931232513621, 'fraction_numerical': 0.03269032543071713, 'mean_word_length': 4.566803278688525, 'pattern_counts': {'":': 0, '<': 108, '<?xml version=': 0, '>': 108, 'https://': 1, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 10, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '24', 'dir': 'jupyter-structured-clean-dedup', 'id': '444371', 'license': 'BSD-2-Clause', 'n_tokens_mistral': 5314, 'n_tokens_neox': 4704, 'n_words': 916, 'path': 'General_Information/.ipynb_checkpoints/glycolysis_itaconate-pwy-zhao2019-checkpoint.ipynb', 'repo_name': 'PNNL-CompBio/bayesian-inference', 'size': '13630'}