It is never easy or cheap to carry out biological research without a lab work, microscopes, etc. Such constraints make challenging to design practical and portable demonstrations for talks or outreach activities. On the other hand, using mathematical models based on physical and chemical principles allow us to investigate these using the dynamical systems machinery and computer simulations, test hypotheses and make predictions about the mechanisms of interest. However, unless the target audience is versed with such concepts, it can be tricky to translate any results using only equations and plots.

In this post, I present examples of systems for which coupling, maths and programming and simple experiments, can be used to present some core concepts to biologists, physicists and general audiences. Specifically how biological oscillators such as neurons and cardiac cells, exhibit similar behaviours to those observed in a simple mechanical system.

Again, for this project the simulations I show here, the numerical codes are written in C and compiled using emscripten. The visualisatios are built using my favorite JS library three.js.

I would like to put under consideration the oscillating system below:

The device is a toy version of the “shishi-odoshi” which means “deer-frigthening” or “boar-frigthening” which according to Wikipedia has those traditional uses. Prettier versions of these can be found in some Japanese gardens.

Beyond beauty and elegance, these devices are interesting because exhibit features shared by some very important biological systems:

The movie below (borrowed from Flavio Fenton’s old website) shows a myocyte being stimulated at a single location. The stimulus location releases calcium which then diffuses through the cell, activating more release locations sites and forming a wave. The stimulus is continuous, yet the behaviour of the myocyte is that of a beating cell.

The diagram above shows several cell compartments connected through pumps and gates (Left panel). The dynamics responsible for the calcium release are shown in the right panel. The Ryanodine receptors (RyR2) can be found in four configurations or states \( (I_1, I_2, O,C) \). The reciprocal of the activation rate \( \tau_r = \frac{1}{k_B} \) corresponds to the recovery from inactivation time.

The translation of the sketched currents in the diagrams into mathematical models can be as complex and detailed as we require (see [1]). In this post however I would like to use the model presented in [2], which is a simplified and elegant version that serves very well to illustrate what I want to discuss here. The diagram of the currents is shown below:

Left: Simplified cell compartment model schematics. Right: The RyRs states Markov model with the four configurations O: open, C: closed, and Inactive states I₁ and I₂.

The equations used in [2] read as follows:

$$\dot{c}_{j} = I_{s}(t)+g_{rel}c_{rel}\frac{k_ac_j^2}{k_{om}+k_ac_j^2}(c_{SR}-c_j) - (c_j-c_o)/\tau_{diff}\quad(1)$$

$$\dot{c}_{rel} = k_{im}(1-c_{rel}) - k_ic_jc_{rel}\quad(2)$$

Eqns. (1) and (2) constitute a minimal model for calcium release. The first term in the rhs of (1) is the external stimulus \(I_{CaL}(t)=I_{s}(t)\), Which is provided by the pacemaker cells in the heart or externally in experiments.

A healthy cell normally releases Calcium periodically, driven by the periodic stimulus and if the recovery time of the ion release channels is fast compared to the stimulus period. An example of this is illustrated below. Where \(T_s=400\ ms\) and \(T_r = 200\ ms \)

If on the other hand, the release is slow due to a malfunction of the calcium release units, the periodic response alternates between a high/low beat to beat calcium release. These calcium alternans are a known mechanism for inducing cardiac arrhythmias such as atrial fibrillation. Below I show results for the case of \( T_s=400\ ms \) and \(T_r = 320\ ms \) .

The complex alternating behaviour of the cardiac cells occurs if the cell possess dysfunctional calcium release units due to a damaged cell or an underlying genetic condition. Mathematically this is encoded in the value of \( k_{im} \).

To investigate this behaviour using the mechanical system, it is easy to build the japanese seesaw with bits and bobs of Lego and K’Nex toys. It however important to take into consideration some of the system important parameters.

To incorporate the stimulus, we need to add mass to the initial mass attached to the right end of the rod. A way to do this is to use the function:

$$m_2=m_{2o} + \alpha t$$

If the angle \(\Phi \geq \Phi^{\dagger} \). Where \(\Phi^{{\dagger}}\) is a parameter specifying the range of action if the stimulus and \(\alpha \) the stimulus strength.

The release condition \( m_2(t') = \gamma m_2(t) \) if \( \Phi = \Phi_2\), where \(\gamma \in [0,1] \) is a parameter stating the fraction of water remaining in the system after the boundary condition at \( \Phi_2 \) is reached. This parameter is the analogous to the parameter \( k_{im}\). It regulates the recovery time.

Other important parameter is the threshold mass \(m_1 \) which needs to be surpassed by \(m_2\) to kickstart the motion.

The release units can be modified in a modular way, and plugged in to the mechanism, just as we do in some gene regulatory networks. The mutants basically consist of release mechanisms which loss water, faster or slower depending on the release position and number of holes in the recipient.

The seesaw in action! This is the WT expected periodic behaviour if the release in almost total and the refill is not very fast, allowing full recovery.

Using one of the mutant release units, we observe the following dynamics:

By replacing the release unit, we observe the 2:1 alternating behaviour. The release unit in this case releases a different amount of water (its the two holed mutant shown above).

Although it quite easy to carry out these type of experiments, It quickly becomes rather complicated to measure a full set of initial conditions and parameter values. For instance to have a cheap and on the fly way to measure the water flow, or to automate the process of water recycling over many tries becomes pretty complicated if we are running the experiments on a near to zero budget. This makes simulation a great tool to use.

Interactive simulator: The tunable parameters are H: Red bar position, H2: Blue bar position, Xo: Stimulus angle range (white line), km: Fraction of water released, Is: current strength. The plots displayed are the mass present at the right end of the seesaw \( \Delta m (t) \) as well as the angle temporal evolution \( \Phi(t) \). Once a solution is computed the PLAY button shows the rod dynamics simulation.

I have used this small system to connect and introduce concepts in genetics, synthetic biology, differential equations, computational physics and recently in physics informed machine learning to audiences of a very varied nature, and thought It would be nice to have it more or less documented in an organised and interactive way accessible for future use.

This is the first out of two entries in the series, the other one linking chemical reactions and electronic circuits.