Derivation of Hill Functions

This notebook follows Example 3.11 in Del Vecchio and Murray, Biomolecular Feedback Systems, and shows how a Hill function can be derived by applying a quasi-steady-state approximation to fast reversible binding reactions.

Model

A transcription factor \(X\) dimerizes before binding a promoter site \(p\). The bound complex \(C\) activates production of the mRNA \(m_Y\), which is translated into protein \(Y\):

\[X + X \xrightleftharpoons[k_2]{k_1} X_2, \qquad X_2 + p \xrightleftharpoons[d]{a} C, \qquad C \xrightarrow{k_f} m_Y + C,\]

\[m_Y \xrightarrow{\kappa} m_Y + Y, \qquad m_Y \xrightarrow{\delta} \emptyset, \qquad Y \xrightarrow{\gamma} \emptyset, \qquad p + C = p_{tot}.\]

Treating \(X(t)\) as an input, the ODE model is

\[\begin{split}\begin{aligned} \frac{dX_2}{dt} &= k_1X^2-k_2X_2-aX_2(p_{tot}-C)+dC, \\ \frac{dC}{dt} &= aX_2(p_{tot}-C)-dC, \\ \frac{dm_Y}{dt} &= k_fC-\delta m_Y, \\ \frac{dY}{dt} &= \kappa m_Y-\gamma Y. \end{aligned}\end{split}\]

from sympy import simplify, symbols

from autoreduce import System, solve_timescale_separation

X, X2, C, mY, Y = symbols("X X2 C m_Y Y")
k1, k2, a, d, kf, delta, kappa, gamma, p_tot = symbols(
    "k1 k2 a d k_f delta kappa gamma p_tot"
)

system = System(
    [X2, C, mY, Y],
    [
        k1 * X**2 - k2 * X2 - a * X2 * (p_tot - C) + d * C,
        a * X2 * (p_tot - C) - d * C,
        kf * C - delta * mY,
        kappa * mY - gamma * Y,
    ],
    params=[X, k1, k2, a, d, kf, delta, kappa, gamma, p_tot],
    params_values=[1.0] * 10,
    x_init=[0.0, 0.0, 0.0, 0.0],
)

Fast binding limit

The binding reactions are much faster than mRNA and protein production and decay. With

\[K_m := \frac{k_2}{k_1}, \qquad K_d := \frac{d}{a}, \qquad c := \frac{k_2}{d}, \qquad \epsilon := \frac{\gamma}{d},\]

the fast equations are written in singular perturbation form by using

\[d=\frac{\gamma}{\epsilon}, \qquad a=\frac{\gamma}{K_d\epsilon}, \qquad k_1=c\frac{\gamma}{K_m\epsilon}, \qquad k_2=c\frac{\gamma}{\epsilon}.\]

Setting \(\epsilon=0\) gives the slow manifold

\[X_2 = \frac{X^2}{K_m}, \qquad C = \frac{p_{tot}X^2/(K_mK_d)}{1 + X^2/(K_mK_d)}.\]

Solve the QSSA reduction

AutoReduce obtains the slow model by retaining \(m_Y\) and \(Y\) and treating \(X_2\) and \(C\) as fast states.

reduced_system, collapsed_system = solve_timescale_separation(
    system,
    slow_states=[mY, Y],
    fast_states=[X2, C],
)

[simplify(expr) for expr in reduced_system.f]

Hill function form

The reduced mRNA dynamics can be written as

\[\frac{dm_Y}{dt} = k_f\frac{p_{tot}X^2/(K_mK_d)}{1+X^2/(K_mK_d)} - \delta m_Y,\]

and the protein equation remains

\[\frac{dY}{dt}=\kappa m_Y-\gamma Y.\]

Letting \(\alpha=k_fp_{tot}\) and \(K=\sqrt{K_mK_d}\) gives the standard Hill expression

\[F(X)=\alpha\frac{(X/K)^2}{1+(X/K)^2}.\]

Km, Kd = symbols("K_m K_d")

production_term = simplify(reduced_system.f[0] + delta * mY)
hill_form = simplify(production_term.subs({k2: Km * k1, d: Kd * a}))

hill_form