HillTau: A fast, compact abstraction for model reduction in biochemical signaling networks

Signaling networks mediate many aspects of cellular function. The conventional, mechanistically motivated approach to modeling such networks is through mass-action chemistry, which maps directly to biological entities and facilitates experimental tests and predictions. However such models are complex, need many parameters, and are computationally costly. Here we introduce the HillTau form for signaling models. HillTau retains the direct mapping to biological observables, but it uses far fewer parameters, and is 100 to over 1000 times faster than ODE-based methods. In the HillTau formalism, the steady-state concentration of signaling molecules is approximated by the Hill equation, and the dynamics by a time-course tau. We demonstrate its use in implementing several biochemical motifs, including association, inhibition, feedforward and feedback inhibition, bistability, oscillations, and a synaptic switch obeying the BCM rule. The major use-cases for HillTau are system abstraction, model reduction, scaffolds for data-driven optimization, and fast approximations to complex cellular signaling.

Reviewer #3: I thank the author for his responses, which largely address my concerns. This manuscript has been strengthened by his edits, especially the addition of a summary of the nature of the HillTau model at the beginning of the results section. In brief, HillTau can be viewed as a phenomenological model that lumps multiple reactions together. Key features are its event-driven nature which allows simulations to make relatively large advances as compared to differential-equation based models, and asymmetric exponential change. Feedback cycles are not directly compatible with the model, but can be approximated by breaking the cycle and using small advances.

Thank you
Although, honestly, I'm still a little confused because the approach uses an internal timestep so it can't skip from event to event. This is a good point, I have clarified in the text around line 178.
For typical use-cases, such as synaptic plasticity models, the event interval is shorter than the time-courses in the model (typically ~1 sec) and hence only a single step is taken. In cases where HillTau inserts additional time-steps for accuracy, it is done behind the scenes of the same event-driven programming interface.
The existing implementation on GitHub is pure Python, but the readme promises a faster C++ version that preserves the existing Python API. This will be invaluable.
A mature C++ version has been on GitHub since July, but it has taken some while to implement a pip-install capability for a range of operating systems. Now this works for Linux platforms and others are in the works. Thank you for pointing this out. I have updated the text to consolidate. Line 196 now leads into the following text, and the duplicate text below has been removed.
MASH runs the reference model through a range of stimuli designed to explore its inputoutput properties, and then uses numerical optimization methods from scipy.optimize to tune parameters so that the HillTau model produces a good fit to the original. I agree. I have inserted the following text to explain.
In brief, FindSim provides a Python-based framework for matching models to experiments. It codifies the experiment design (e.g., time-series, dose-response, bar-chart) and experimental results into a single machine-readable file. FindSim runs the experiment on the model and returns a numerical score for goodness of fit. The model may be defined in SBML (run using MOOSE) or using HillTau. Thus FindSim can be used as the scoring function for optimizing model fit to experiments using a variety of optimization methods available in scipy.optimize. lines 457-458: not 100% clear what you're considering state variables here; usually these would be things that have associated rates of change, and then you very much would want to know the initial values. I'm thinking of things like gating variables in Hodgkin-Huxley. There's no initial concentration but there's definitely an initial value that matters. This is a good point. I now address it in the text: Here we consider state variables to be those which are computed, as opposed to defined using initial conditions. In ODE models we estimate this by counting the number of rate terms plus the number of molecular species with a non-zero initial value. In HillTau models we count the rate terms and the species which are not reaction outputs. This yields the approximate scaling terms below.
Line 469: mysterious space in the middle of C++.

Fixed.
Line 914: I think "HillTaul" is a typo, otherwise I missed something major.

Fixed
In addition to the above changes, I have updated the benchmark plots in Figure 6 to reflect further improvements to the automated timestep algorithm in HillTau. The regression fit scores in the figure 6 legend and text have been updated to match.