Predict gene knockout strategies

In cameo we have two ways of predicting gene knockout targets: using evolutionary algorithms (OptGene) or linear programming (OptKnock)

If you’re running this notebook on try.cameo.bio, things might run very slow due to our inability to provide access to the proprietary CPLEX solver on a public webserver. Furthermore, Jupyter kernels might crash and restart due to memory limitations on the server.

from cameo import models
No module named ggplot
/Users/niso/Dev/cameo/cameo/visualization/plotting/__init__.py:53 UserWarning: Cannot import any plotting library. Please install one of 'plotly', 'bokeh' or 'ggplot' if you want to use any plotting function.
from cameo.visualization.plotting import plotter
plotter._engine
<cameo.visualization.plotting.abstract.AbstractPlotter at 0x108debc10>
model = models.bigg.iJO1366
wt_solution = model.optimize()
growth = wt_solution.fluxes["BIOMASS_Ec_iJO1366_core_53p95M"]
acetate_production = wt_solution.fluxes["EX_ac_e"]
from cameo import phenotypic_phase_plane
p = phenotypic_phase_plane(model, variables=['BIOMASS_Ec_iJO1366_core_53p95M'], objective='EX_ac_e')
# p.plot(points=[(growth, acetate_production)])
p.data_frame
BIOMASS_Ec_iJO1366_core_53p95M objective_lower_bound objective_upper_bound c_yield_lower_bound c_yield_upper_bound mass_yield_lower_bound mass_yield_upper_bound
0 0.000000 0.0 2.909347e+01 0.0 9.697822e-01 0.0 9.535050e-01
1 0.051704 0.0 2.770228e+01 0.0 9.234092e-01 0.0 9.079103e-01
2 0.103408 0.0 2.631109e+01 0.0 8.770362e-01 0.0 8.623157e-01
3 0.155111 0.0 2.491990e+01 0.0 8.306632e-01 0.0 8.167210e-01
4 0.206815 0.0 2.352871e+01 0.0 7.842902e-01 0.0 7.711263e-01
5 0.258519 0.0 2.213752e+01 0.0 7.379172e-01 0.0 7.255317e-01
6 0.310223 0.0 2.074633e+01 0.0 6.915442e-01 0.0 6.799370e-01
7 0.361926 0.0 1.935220e+01 0.0 6.450733e-01 0.0 6.342461e-01
8 0.413630 0.0 1.795528e+01 0.0 5.985094e-01 0.0 5.884637e-01
9 0.465334 0.0 1.655358e+01 0.0 5.517861e-01 0.0 5.425247e-01
10 0.517038 0.0 1.509136e+01 0.0 5.030452e-01 0.0 4.946018e-01
11 0.568742 0.0 1.362913e+01 0.0 4.543042e-01 0.0 4.466790e-01
12 0.620445 0.0 1.216690e+01 0.0 4.055633e-01 0.0 3.987561e-01
13 0.672149 0.0 1.069047e+01 0.0 3.563490e-01 0.0 3.503679e-01
14 0.723853 0.0 9.191568e+00 0.0 3.063856e-01 0.0 3.012431e-01
15 0.775557 0.0 7.493955e+00 0.0 2.497985e-01 0.0 2.456058e-01
16 0.827260 0.0 5.753441e+00 0.0 1.917814e-01 0.0 1.885624e-01
17 0.878964 0.0 3.835628e+00 0.0 1.278543e-01 0.0 1.257083e-01
18 0.930668 0.0 1.917814e+00 0.0 6.392713e-02 0.0 6.285414e-02
19 0.982372 0.0 3.957945e-13 0.0 1.319315e-14 0.0 1.297171e-14
import altair
c = altair.Chart(p.data_frame)
c.mark_line().encode(x='BIOMASS_Ec_iJO1366_core_53p95M', y='objective_upper_bound')
05-predict-gene-knockout-strategies_files/05-predict-gene-knockout-strategies_12_2.png
import altair as alt

# load built-in dataset as a pandas DataFrame
cars = alt.load_dataset('cars')

# Uncomment for rendering in JupyterLab & nteract
# alt.enable_mime_rendering()

alt.Chart(cars).mark_circle().encode(
    x='Horsepower',
    y='Miles_per_Gallon',
    color='Origin',
)
<altair.VegaLite object>

OptGene

OptGene is an approach to search for gene or reaction knockouts that relies on evolutionary algorithms[1]. The following image from authors summarizes the OptGene workflow.

Every iteration we keep the best 50 individuals so we can generate a library of targets.

from cameo.strain_design.heuristic.evolutionary_based import OptGene
optgene = OptGene(model)
result = optgene.run(target=model.reactions.EX_ac_e,
                     biomass=model.reactions.BIOMASS_Ec_iJO1366_core_53p95M,
                     substrate=model.metabolites.glc__D_e,
                     max_evaluations=5000,
                     plot=False)
Starting optimization at Wed, 25 Oct 2017 11:18:50
/Users/niso/.local/share/virtualenvs/cameo2.7/lib/python2.7/site-packages/bokeh/models/sources.py:137: BokehUserWarning:

ColumnDataSource's columns must be of the same length. Current lengths: ('x', 1), ('y', 0)

Failed to display Jupyter Widget of type HBox.

If you're reading this message in Jupyter Notebook or JupyterLab, it may mean that the widgets JavaScript is still loading. If this message persists, it likely means that the widgets JavaScript library is either not installed or not enabled. See the Jupyter Widgets Documentation for setup instructions.

If you're reading this message in another notebook frontend (for example, a static rendering on GitHub or NBViewer), it may mean that your frontend doesn't currently support widgets.

/Users/niso/.local/share/virtualenvs/cameo2.7/lib/python2.7/site-packages/bokeh/models/sources.py:137: BokehUserWarning:

ColumnDataSource's columns must be of the same length. Current lengths: ('x', 2), ('y', 1)

/Users/niso/.local/share/virtualenvs/cameo2.7/lib/python2.7/site-packages/bokeh/models/sources.py:137: BokehUserWarning:

ColumnDataSource's columns must be of the same length. Current lengths: ('x', 3), ('y', 2)

/Users/niso/.local/share/virtualenvs/cameo2.7/lib/python2.7/site-packages/bokeh/models/sources.py:137: BokehUserWarning:

ColumnDataSource's columns must be of the same length. Current lengths: ('x', 4), ('y', 3)

/Users/niso/.local/share/virtualenvs/cameo2.7/lib/python2.7/site-packages/bokeh/models/sources.py:137: BokehUserWarning:

ColumnDataSource's columns must be of the same length. Current lengths: ('x', 5), ('y', 4)

/Users/niso/.local/share/virtualenvs/cameo2.7/lib/python2.7/site-packages/bokeh/models/sources.py:137: BokehUserWarning:

ColumnDataSource's columns must be of the same length. Current lengths: ('x', 6), ('y', 5)

/Users/niso/.local/share/virtualenvs/cameo2.7/lib/python2.7/site-packages/bokeh/models/sources.py:137: BokehUserWarning:

ColumnDataSource's columns must be of the same length. Current lengths: ('x', 7), ('y', 6)

/Users/niso/.local/share/virtualenvs/cameo2.7/lib/python2.7/site-packages/bokeh/models/sources.py:137: BokehUserWarning:

ColumnDataSource's columns must be of the same length. Current lengths: ('x', 8), ('y', 7)

/Users/niso/.local/share/virtualenvs/cameo2.7/lib/python2.7/site-packages/bokeh/models/sources.py:137: BokehUserWarning:

ColumnDataSource's columns must be of the same length. Current lengths: ('x', 9), ('y', 8)

/Users/niso/.local/share/virtualenvs/cameo2.7/lib/python2.7/site-packages/bokeh/models/sources.py:137: BokehUserWarning:

ColumnDataSource's columns must be of the same length. Current lengths: ('x', 10), ('y', 9)

/Users/niso/.local/share/virtualenvs/cameo2.7/lib/python2.7/site-packages/bokeh/models/sources.py:137: BokehUserWarning:

ColumnDataSource's columns must be of the same length. Current lengths: ('x', 11), ('y', 10)

/Users/niso/.local/share/virtualenvs/cameo2.7/lib/python2.7/site-packages/bokeh/models/sources.py:137: BokehUserWarning:

ColumnDataSource's columns must be of the same length. Current lengths: ('x', 12), ('y', 11)

/Users/niso/.local/share/virtualenvs/cameo2.7/lib/python2.7/site-packages/bokeh/models/sources.py:137: BokehUserWarning:

ColumnDataSource's columns must be of the same length. Current lengths: ('x', 13), ('y', 12)

/Users/niso/.local/share/virtualenvs/cameo2.7/lib/python2.7/site-packages/bokeh/models/sources.py:137: BokehUserWarning:

ColumnDataSource's columns must be of the same length. Current lengths: ('x', 14), ('y', 13)

/Users/niso/.local/share/virtualenvs/cameo2.7/lib/python2.7/site-packages/bokeh/models/sources.py:137: BokehUserWarning:

ColumnDataSource's columns must be of the same length. Current lengths: ('x', 15), ('y', 14)

/Users/niso/.local/share/virtualenvs/cameo2.7/lib/python2.7/site-packages/bokeh/models/sources.py:137: BokehUserWarning:

ColumnDataSource's columns must be of the same length. Current lengths: ('x', 16), ('y', 15)

/Users/niso/.local/share/virtualenvs/cameo2.7/lib/python2.7/site-packages/bokeh/models/sources.py:137: BokehUserWarning:

ColumnDataSource's columns must be of the same length. Current lengths: ('x', 17), ('y', 16)

/Users/niso/.local/share/virtualenvs/cameo2.7/lib/python2.7/site-packages/bokeh/models/sources.py:137: BokehUserWarning:

ColumnDataSource's columns must be of the same length. Current lengths: ('x', 18), ('y', 17)

/Users/niso/.local/share/virtualenvs/cameo2.7/lib/python2.7/site-packages/bokeh/models/sources.py:137: BokehUserWarning:

ColumnDataSource's columns must be of the same length. Current lengths: ('x', 19), ('y', 18)

/Users/niso/.local/share/virtualenvs/cameo2.7/lib/python2.7/site-packages/bokeh/models/sources.py:137: BokehUserWarning:

ColumnDataSource's columns must be of the same length. Current lengths: ('x', 20), ('y', 19)

/Users/niso/.local/share/virtualenvs/cameo2.7/lib/python2.7/site-packages/bokeh/models/sources.py:137: BokehUserWarning:

ColumnDataSource's columns must be of the same length. Current lengths: ('x', 21), ('y', 20)

/Users/niso/.local/share/virtualenvs/cameo2.7/lib/python2.7/site-packages/bokeh/models/sources.py:137: BokehUserWarning:

ColumnDataSource's columns must be of the same length. Current lengths: ('x', 22), ('y', 21)

/Users/niso/.local/share/virtualenvs/cameo2.7/lib/python2.7/site-packages/bokeh/models/sources.py:137: BokehUserWarning:

ColumnDataSource's columns must be of the same length. Current lengths: ('x', 23), ('y', 22)

/Users/niso/.local/share/virtualenvs/cameo2.7/lib/python2.7/site-packages/bokeh/models/sources.py:137: BokehUserWarning:

ColumnDataSource's columns must be of the same length. Current lengths: ('x', 24), ('y', 23)

/Users/niso/.local/share/virtualenvs/cameo2.7/lib/python2.7/site-packages/bokeh/models/sources.py:137: BokehUserWarning:

ColumnDataSource's columns must be of the same length. Current lengths: ('x', 25), ('y', 24)
Finished after 00:14:20
result

Failed to display Jupyter Widget of type HBox.

If you're reading this message in Jupyter Notebook or JupyterLab, it may mean that the widgets JavaScript is still loading. If this message persists, it likely means that the widgets JavaScript library is either not installed or not enabled. See the Jupyter Widgets Documentation for setup instructions.

If you're reading this message in another notebook frontend (for example, a static rendering on GitHub or NBViewer), it may mean that your frontend doesn't currently support widgets.

OptGene Result

  • Simulation: fba
  • Objective Function: $$bpcy = \frac{(BIOMASS\_Ec\_iJO1366\_core\_53p95M * EX\_ac\_e)}{EX\_glc\_\_D\_e}$$
reactions genes size fva_min fva_max target_flux biomass_flux yield fitness
0 (ATPS4rpp,) ((b3736,), (b3733,), (b3732,), (b3731,), (b373... 1 0.0 14.187817 13.942930 0.402477 1.394293 0.561171
1 (ATPS4rpp, PGCD) ((b3732, b2913),) 2 0.0 14.976296 14.801391 0.388364 1.480139 0.574833
2 (ATPS4rpp, PSP_L) ((b3731, b4388),) 2 0.0 14.976296 14.801391 0.388364 1.480139 0.574833
result.plot(0)
result.display_on_map(0, "iJO1366.Central metabolism")

OptKnock

OptKnock uses a bi-level mixed integer linear programming approach to identify reaction knockouts[2]:

\[\begin{split}\begin{matrix} maximize & \mathit{v_{chemical}} & & (\mathbf{OptKnock}) \\ \mathit{y_j} & & & \\ subject~to & maximize & \mathit{v_{biomass}} & (\mathbf{Primal}) \\ & \mathit{v_j} & & & & \\ \end{matrix}\\ \begin{bmatrix} subject~to & \sum_{j=1}^{M}S_{ij}v_{j} = 0,\\ & v_{carbon\_uptake} = v_{carbon~target}\\ & v_{apt} \ge v_{apt\_main}\\ & v_{biomass} \ge v_{target\_biomass}\\ & v_{j}^{min} \cdot y_j \le v_j \le v_{j}^{max} \cdot y_j, \forall j \in \boldsymbol{M} \\ \end{bmatrix}\\ \begin{align} & y_j = {0, 1}, & & \forall j \in \boldsymbol{M} & \\ & \sum_{j \in M} (1 - y_j) \le K& & & \\ \end{align}\end{split}\]
from cameo.strain_design.deterministic.linear_programming import OptKnock
optknock = OptKnock(model, fraction_of_optimum=0.1)
/Users/niso/Dev/cameo/cameo/strain_design/deterministic/linear_programming.py:104: UserWarning:

You are trying to run OptKnock with glpk_interface. This might not end well.
---------------------------------------------------------------------------

KeyboardInterrupt                         Traceback (most recent call last)

<ipython-input-14-f9841f906563> in <module>()
----> 1 optknock = OptKnock(model, fraction_of_optimum=0.1)


/Users/niso/Dev/cameo/cameo/strain_design/deterministic/linear_programming.pyc in __init__(self, model, exclude_reactions, remove_blocked, fraction_of_optimum, exclude_non_gene_reactions, use_nullspace_simplification, *args, **kwargs)
    107             fix_objective_as_constraint(self._model, fraction=fraction_of_optimum)
    108         if remove_blocked:
--> 109             self._remove_blocked_reactions()
    110         if not exclude_reactions:
    111             exclude_reactions = []


/Users/niso/Dev/cameo/cameo/strain_design/deterministic/linear_programming.pyc in _remove_blocked_reactions(self)
    116
    117     def _remove_blocked_reactions(self):
--> 118         fva_res = flux_variability_analysis(self._model, fraction_of_optimum=0)
    119         blocked = [
    120             self._model.reactions.get_by_id(reaction) for reaction, row in fva_res.data_frame.iterrows()


/Users/niso/Dev/cameo/cameo/flux_analysis/analysis.pyc in flux_variability_analysis(model, reactions, fraction_of_optimum, pfba_factor, remove_cycles, view)
    221         else:
    222             func_obj = _FvaFunctionObject(model, _flux_variability_analysis)
--> 223         chunky_results = view.map(func_obj, reaction_chunks)
    224         solution = pandas.concat(chunky_results)
    225


/Users/niso/Dev/cameo/cameo/parallel.pyc in map(self, *args, **kwargs)
    231 class SequentialView(object):
    232     def map(self, *args, **kwargs):
--> 233         return list(map(*args, **kwargs))
    234
    235     def apply(self, func, *args, **kwargs):


/Users/niso/Dev/cameo/cameo/flux_analysis/analysis.pyc in __call__(self, reactions)
    343
    344     def __call__(self, reactions):
--> 345         return self.fva(self.model, reactions)
    346
    347


/Users/niso/Dev/cameo/cameo/flux_analysis/analysis.pyc in _flux_variability_analysis(model, reactions)
    362             model.solver.objective.set_linear_coefficients({reaction.forward_variable: 1.,
    363                                                             reaction.reverse_variable: -1.})
--> 364             model.solver.optimize()
    365             if model.solver.status == OPTIMAL:
    366                 fva_sol[reaction.id]['lower_bound'] = model.objective.value


/Users/niso/.local/share/virtualenvs/cameo2.7/lib/python2.7/site-packages/optlang/interface.pyc in optimize(self)
   1441         """
   1442         self.update()
-> 1443         status = self._optimize()
   1444         if status != OPTIMAL and self.configuration.presolve == "auto":
   1445             self.configuration.presolve = True


/Users/niso/.local/share/virtualenvs/cameo2.7/lib/python2.7/site-packages/optlang/glpk_interface.pyc in _optimize(self)
    670
    671     def _optimize(self):
--> 672         status = self._run_glp_simplex()
    673
    674         # Sometimes GLPK gets itself stuck with an invalid basis. This will help it get rid of it.


/Users/niso/.local/share/virtualenvs/cameo2.7/lib/python2.7/site-packages/optlang/glpk_interface.pyc in _run_glp_simplex(self)
    644
    645     def _run_glp_simplex(self):
--> 646         return_value = glp_simplex(self.problem, self.configuration._smcp)
    647         glpk_status = glp_get_status(self.problem)
    648         if return_value == 0:


KeyboardInterrupt:

Running multiple knockouts with OptKnock can take a few hours or days…

result = optknock.run(max_knockouts=1, target="EX_ac_e", biomass="BIOMASS_Ec_iJO1366_core_53p95M")
result
result.plot(0)
result.display_on_map(0, "iJO1366.Central metabolism")

References

[1]Patil, K. R., Rocha, I., Förster, J., & Nielsen, J. (2005). Evolutionary programming as a platform for in silico metabolic engineering. BMC Bioinformatics, 6, 308. doi:10.1186/1471-2105-6-308

[2]Burgard, A.P., Pharkya, P., Maranas, C.D. (2003), “OptKnock: A Bilevel Programming Framework for Identifying Gene Knockout Strategies for Microbial Strain Optimization,” Biotechnology and Bioengineering, 84(6), 647-657.