Python · chapter 7 of 13

Time-varying stimuli

Real experiments switch genes on, wash inducers out, ramp a dose, or drive a circuit periodically. grn.stimulus describes such an input and grn.simulate (or grn.stimulate, grn.steady_state, grn.knockout_screen) applies it.

How a stimulus enters the model

A stimulus acts on the production (transcription) rate of one or more genes through a waveform u(t):

`mode`	Production of gene i	Interpretation
`"scale"` (default)	`u(t) * production_i(x)`	An inducible promoter: `u = 0` silences the gene, `u = 1` is the unstimulated rate, `u = 3` triples it
`"add"`	`production_i(x) + u(t)`	An external source of transcript in production units

Scale stimuli on the same gene multiply, add stimuli sum. Regulation by the network still applies: a stimulated gene with repressors bound still produces little.

Waveforms

s = grn.stimulus(gene, type, ...) where gene is an index, a label, or a list of either.

`type`	`u(t)`	Arguments
`"step"`	`amplitude` for `start <= t < stop`, `baseline` elsewhere	`start` (0), `stop` (inf), `rise` (0)
`"ramp"`	`baseline` before `start`, linear to `amplitude` at `stop`, `amplitude` after	`start`, `stop` (required)
`"sine"`	`baseline + (amplitude - baseline) (1 - cos(2 pi (t - start) / period)) / 2` for `start <= t < stop`	`period` (10), `start`, `stop`
`"custom"`	`function(t)`	`function` (required), `breaks`

Defaults: baseline=1 and amplitude=2 for scale mode, baseline=0 and amplitude=1 for add mode. A step with a finite stop is a pulse. rise > 0 replaces the jumps of a step by logistic ramps with roughly that 10–90 % rise time, for experiments where the inducer diffuses in.

on     = grn.stimulus("X", "step", baseline=0, amplitude=1, start=5)
pulse  = grn.stimulus(0, "step", amplitude=3, start=5, stop=10)
dose   = grn.stimulus(0, "ramp", baseline=0, amplitude=3, start=0, stop=40)
drive  = grn.stimulus("G2", "sine", amplitude=3, period=12)
custom = grn.stimulus(0, "custom", function=lambda t: 1 + 0.5 * (t > 3) + 0.5 * (t > 6), breaks=[3, 6])

The built-in waveforms. A step with a finite stop is a pulse, and rise turns its edges into logistic ramps.

Timing and the solver

Stimulus times are on the clock of the recorded solution: t = 0 is the first recorded sample and a transient runs at negative times. Every step edge, ramp corner, sine window edge, and custom breaks time is passed to the ODE integrator, which restarts there. A gene switched on at t = 4 is therefore exactly zero before t = 4, not leaked across a solver step. The Euler–Maruyama integrator evaluates the stimulus at every step.

Applying stimuli

sol = grn.simulate(model, stimulus=on, duration=30)             # one
sol = grn.simulate(model, stimulus=[on, drive], duration=60)    # several
model2 = grn.stimulate(model, pulse)      # attach to a copy; model2.f(t, x) includes it
model2 = grn.stimulate(model2, None)      # remove all stimuli

The realised input is recorded in sol.u (the multiplicative factor, (samples, n), all ones without stimuli) and sol.u_add (the additive term). grn.plot_expression adds a panel showing the input when the solution carries a stimulus.

np.column_stack([sol.t[47:53], sol.u[47:53]])      # the input around the edge at t = 5

Output

[[4.7 0.  1.  1. ]
 [4.8 0.  1.  1. ]
 [4.9 0.  1.  1. ]
 [5.  1.  1.  1. ]
 [5.1 1.  1.  1. ]
 [5.2 1.  1.  1. ]]

Steady states under a stimulus

Steady states need an autonomous system. grn.steady_state and grn.knockout_screen freeze the stimulus at its value at time (default 0):

off = grn.steady_state(model, stimulus=on, time=0)    # before the step
xon = grn.steady_state(model, stimulus=on, time=10)   # after it
ko  = grn.knockout_screen(model, stimulus=on, time=10)

This is how to build dose-response steady-state data: a step of each amplitude, time anywhere after start.

model = grn.build_model(grn.motif("sim", size=5), seed=1)
doses = np.linspace(0, 3, 13)
X = np.array([grn.steady_state(model, stimulus=grn.stimulus("TF", "step", baseline=d, amplitude=d))
              for d in doses])

A single-input module driven from 0 to 3 times its normal rate. Left: steady-state dose response of the four targets, each with its own threshold and Hill coefficient. Right: the same information read dynamically from a slow ramp.

Example: the incoherent feed-forward loop pulse

model = grn.build_model("ffl_incoherent", logic="multiplicative", seed=1,
                        params=dict(gamma=[5, 1, 3], vmax=5, basal=[5, 0, 0]))
on  = grn.stimulus("X", "step", baseline=0, amplitude=1, start=5)
sol = grn.simulate(model, stimulus=on, initial_condition=np.zeros(3), duration=30)
grn.plot_expression(sol)     # Z rises when X switches on, then Y catches up and represses it

The incoherent feed-forward loop with X switched on at t = 5. The top panel shows the input; Z pulses because X activates it directly before Y, which X also activates, arrives to repress it.

← Perturbations, steady states, and knockout screens Figures →