Here’s the puzzle a metabolic engineer faces. You have E. coli. You want it to make, say, lots of succinate. You sketch out the metabolism on a whiteboard — glucose comes in, gets chopped through glycolysis, enters the TCA cycle, and somewhere in there is the reaction you care about. Now: if you knock out gene A and overexpress gene B, what will the succinate yield be?
The textbook answer is: simulate the kinetics. Every enzyme has a rate law (Michaelis–Menten, or worse), and if you know all the $V_{\max}$, all the $K_m$, all the enzyme concentrations, you can write down a giant system of ODEs and integrate it.
The trouble is you don’t know any of those numbers. Not really. The published values are wildly inconsistent, vary with conditions, and don’t measure the same thing in the cell as in the test tube. For E. coli you’ve got maybe a few hundred enzymes with okay-ish kinetic data, out of a couple thousand reactions. For non-model organisms, you’ve got essentially nothing.
So how can we predict anything?
Forget rate constants for a moment. The cell still has to obey mass balance. If a metabolite is being produced and consumed inside the cell, and the cell is in steady state (not piling up the metabolite, not running out), then production has to equal consumption. That’s it. That’s the law.
Write each reaction as a vector of stoichiometric coefficients. Glycolysis’s first reaction, hexokinase, is:
\[\text{glucose} + \text{ATP} \to \text{glucose-6-P} + \text{ADP}\]So in the column for hexokinase, you’d have $-1$ for glucose, $-1$ for ATP, $+1$ for G6P, $+1$ for ADP, and $0$ for everything else.
Stack all those columns side by side, one per reaction. You get a matrix $S$. Rows are metabolites, columns are reactions. For E. coli, $S$ is roughly $2000 \times 2500$. Sparse, mostly zeros.
Now let $v$ be a vector of fluxes — how fast each reaction runs, in mmol per gram dry weight per hour. The change in each metabolite’s concentration is:
\[\frac{dc}{dt} = S v\]At steady state, $dc/dt = 0$. So:
\[S v = 0\]That’s the whole governing equation. It is a linear system. And it has way more unknowns (reactions) than equations (metabolites), which means there’s an entire subspace of valid flux distributions — a hyperplane in $\mathbb{R}^{2500}$ of states the cell could be in without violating mass balance.
Two more constraints, both physically obvious:
Add these inequality constraints and the feasible flux space becomes a convex polytope. Not a hyperplane anymore — a finite, bounded region in flux-space. Every point in it is a flux distribution that the cell could in principle adopt without breaking mass balance or running enzymes backwards.
The polytope can still be huge. Millions of feasible flux distributions for a real genome-scale model. Which one does the cell actually pick?
This is where flux balance analysis (FBA) makes its bet. The bet is: evolution has shaped E. coli on rich glucose medium to maximize biomass production rate. So the cell picks the point in the polytope that maximizes a biomass objective: a pseudo-reaction that drains the precursors (amino acids, nucleotides, lipids) in the right proportions to build one gram of cell.
Now your problem is:
\[\max_{v} \; c^T v \quad \text{subject to} \quad S v = 0, \; v^{\min} \leq v \leq v^{\max}\]where $c$ is a vector picking out the biomass reaction (or whatever else you want — succinate secretion, say). Linear objective, linear constraints. That is a linear program. We have solved linear programming. cobrapy does it in milliseconds on a laptop.
And — this is the part that should feel a little miraculous — for E. coli on glucose, this predicts growth rate within ~10% of experiment, and predicts most byproduct secretion patterns correctly, without using a single rate constant. The mass balance and the optimization assumption do all the work.
This is where it earns its keep in engineering:
The optimization hypothesis is the soft part. E. coli on glucose evolved to grow fast, sure. E. coli engineered to make succinate did not evolve for that. The cell is happy to lose your succinate flux if it can grow faster. So “optimize biomass” predicts what the wild-type does, which is exactly not what you want.
The standard fixes:
The hard part is that FBA tells you nothing about rates over time, enzyme concentrations, or regulation. A cell can satisfy mass balance with two completely different enzyme expression patterns; FBA cannot distinguish them. If the answer to your question is “this enzyme is rate-limiting because it’s barely expressed,” FBA shrugs. You need 13C metabolic flux analysis (which measures actual fluxes via isotope labeling) or proteomics or kinetic modeling to fill in.
The thing I want you to feel is how much you get for free. A genome-scale stoichiometry matrix is just bookkeeping — you read off the reactions from the genome annotation and assemble $S$ row by row. No experiments needed. Then linear algebra plus a single optimization assumption gives you growth rate, yield ceilings, knockout predictions, and media sensitivity. For a system that has billions of years of evolutionary fine-tuning baked into it.
That’s the move. When you can’t model the dynamics, model the constraints, and bet that evolution has pushed the system to a corner of the feasible polytope. It’s a paradigm worth stealing for other domains — economics has been doing essentially the same thing with input-output models since Leontief in the 1930s, and the resemblance is not coincidental.
The natural follow-up is 13C metabolic flux analysis — how you actually measure in-vivo fluxes (versus predicting them with FBA), by feeding cells labeled glucose and tracing where the labels end up. Pair it with Antoniewicz’s reviews. After that, thermodynamic FBA (Henry et al.) adds Gibbs-free-energy constraints to rule out flux distributions that are stoichiometrically valid but thermodynamically backwards — a substantial tightening that pure FBA misses.