2022-10-17, Klein, Julia, Petrov, Tatjana

Biological groups exhibit fascinating collective dynamics without centralised control, through only local interactions between individuals. Desirable group behaviours are typically linked to a certain fitness function, which the group robustly performs under different perturbations in, for instance, group structure, group size, noise, or environmental factors. Deriving this fitness function is an important step towards understanding the collective response, yet it easily becomes non-trivial in the context of complex collective dynamics. In particular, understanding the social feedback - how the collective behaviour adapts to changes in the group size - requires dealing with complex models and limited experimental data. In this work, we assume that the collective response is experimentally observed for a chosen, finite set of group sizes. Based on such data, we propose a framework which allows to: (i) predict the collective response for any given group size, and (ii) automatically propose a fitness function. We use Smoothed Model Checking, an approach based on Gaussian Process Classification, to develop a methodology that is scalable, flexible, and data-efficient; We specify the fitness function as a template temporal logic formula with unknown parameters, and we automatically infer the missing quantities from data. We evaluate the framework over a case study of a collective stinging defence mechanism in honeybee colonies.

2021-12, Repin, Denis, Petrov, Tatjana

Predicting stochastic cellular dynamics as emerging from the mechanistic models of molecular interactions is a long-standing challenge in systems biology: low-level chemical reaction network (CRN) models give rise to a highly-dimensional continuous-time Markov chain (CTMC) which is computationally demanding and often prohibitive to analyse in practice. A recently proposed abstraction method uses deep learning to replace this CTMC with a discrete-time continuous-space process, by training a mixture density deep neural network with traces sampled at regular time intervals (which can be obtained either by simulating a given CRN or as time-series data from experiment). The major advantage of such abstraction is that it produces a computational model that is dramatically cheaper to execute, while it preserves the statistical features of the training data. In general, the abstraction accuracy improves with the amount of training data. However, depending on the CRN, the overall quality of the method – the efficiency gain and abstraction accuracy – will also depend on the choice of neural network architecture given by hyper-parameters such as the layer types and connections between them. As a consequence, in practice, the modeller has to take care of finding the suitable architecture manually, for each given CRN, through a tedious and time-consuming trial-and-error cycle. In this paper, we propose to further automatise deep abstractions for stochastic CRNs, through learning the neural network architecture along with learning the transition kernel of the abstract process. Automated search of the architecture makes the method applicable directly to any given CRN, which is time-saving for deep learning experts and crucial for non-specialists. We implement the method and demonstrate its performance on a number of representative CRNs with multi-modal emergent phenotypes. Moreover, we showcase that deep abstractions can be used for efficient multi-scale simulations, which are otherwise computationally intractable. To this end, we define a scenario where multiple CRN instances interact across a spatial grid via shared species. Finally, we discuss the limitations and challenges arising when using deep abstractions.

2021, Hajnal, Matej, Šafránek, David, Petrov, Tatjana

We present a tool for inferring the parameters of a Discrete-time Markov chain (DTMC) with respect to properties written in probabilistic temporal logic (PCTL) informed by data observations. The tool combines, in a modular and user-friendly way, the existing methods and tools for parameter synthesis of DTMCs. On top of this, the tool implements several hybrid methods for the exploration of the parameter space based on utilising the intermediate results of parametric model checking – the symbolic representation of properties’ satisfaction in the form of rational functions. These methods are combined to support three different parameter exploration methods: (i) optimisation, (ii) parameter synthesis, (iii) Bayesian parameter inference. Each of the available methods makes a different trade-off between scalability and inference quality, which can be chosen by the user depending on the application context. In this paper, we present the implementation, the main features of the tool, and we evaluate its performance on several benchmarks.

2020, Repin, Denis, Phung, Nhat-Huy, Petrov, Tatjana

We present a toolbox for stochastic simulations with CRN models and their (automated) deep abstractions: a mixture density deep neural network trained on time-series data produced by the CRN. The optimal neural network architecture is learnt along with learning the transition kernel of the abstract process. Automated search of the architecture makes the method applicable directly to any given CRN, which is time-saving for deep learning experts and crucial for non-specialists. The tool was primarily designed to efficiently reproduce simulation traces of given complex stochastic reaction networks arising in systems biology research, possibly with multi-modal emergent phenotypes. It is at the same time applicable to any other application domain, where time-series measurements of a Markovian stochastic process are available by experiment or synthesised with simulation (e.g. are obtained from a rule-based description of the CRN).

2022, Petrov, Tatjana, Hajnal, Matej, Klein, Julia, Šafránek, David, Nouvian, Morgane

Honeybees protect their colony against vertebrates by mass stinging and they coordinate their actions during this crucial event thanks to an alarm pheromone carried directly on the stinger, which is therefore released upon stinging. The pheromone then recruits nearby bees so that more and more bees participate in the defence. However, a quantitative understanding of how an individual bee adapts its stinging response during the course of an attack is still a challenge: Typically, only the group behaviour is effectively measurable in experiment; Further, linking the observed group behaviour with individual responses requires a probabilistic model enumerating a combinatorial number of possible group contexts during the defence; Finally, extracting the individual characteristics from group observations requires novel methods for parameter inference.
We first experimentally observed the behaviour of groups of bees confronted with a fake predator inside an arena and quantified their defensive reaction by counting the number of stingers embedded in the dummy at the end of a trial. We propose a biologically plausible model of this phenomenon, which transparently links the choice of each individual bee to sting or not, to its group context at the time of the decision. Then, we propose an efficient method for inferring the parameters of the model from the experimental data. Finally, we use this methodology to investigate the effect of group size on stinging initiation and alarm pheromone recruitment.
Our findings shed light on how the social context influences stinging behaviour, by quantifying how the alarm pheromone concentration level affects the decision of each bee to sting or not in a given group size. We show that recruitment is curbed as group size grows, thus suggesting that the presence of nestmates is integrated as a negative cue by individual bees. Moreover, the unique integration of exact and statistical methods provides a quantitative characterisation of uncertainty associated to each of the inferred parameters.

2021-11, Petrov, Tatjana, Igler, Claudia, Sezgin, Ali, Henzinger, Thomas A., Guet, Calin C.

Gene expression is regulated by the set of transcription factors (TFs) that bind to the promoter. The ensuing regulating function is often represented as a combinational logic circuit, where output (gene expression) is determined by current input values (promoter bound TFs) only. However, the simultaneous arrival of TFs is a strong assumption, since transcription and translation of genes introduce intrinsic time delays and there is no global synchronization among the arrival times of different molecular species at their targets. We present an experimentally implementable genetic circuit with two inputs and one output, which in the presence of small delays in input arrival, exhibits qualitatively distinct population-level phenotypes, over timescales that are longer than typical cell doubling times. From a dynamical systems point of view, these phenotypes represent long-lived transients: although they converge to the same value eventually, they do so after a very long time span. The key feature of this toy model genetic circuit is that, despite having only two inputs and one output, it is regulated by twenty-three distinct DNA-TF configurations, two of which are more stable than others (DNA looped states), one promoting and another blocking the expression of the output gene. Small delays in input arrival time result in a majority of cells in the population quickly reaching the stable state associated with the first input, while exiting of this stable state occurs at a slow timescale. In order to mechanistically model the behavior of this genetic circuit, we used a rule-based modelling language, and implemented a grid-search to find parameter combinations giving rise to long-lived transients. Our analysis shows that in the absence of feedback, there exist path-dependent gene regulatory mechanisms based on the long timescale of transients. The behavior of this toy model circuit suggests that gene regulatory networks can exploit event timing to create phenotypes, and it opens the possibility that they could use event timing to memorize events, without regulatory feedback. The model reveals the importance of (i) mechanistically modelling the transitions between the different DNA-TF states, and (ii) employing transient analysis thereof.

2020-09-29, Bokes, Pavol, Klein, Julia, Petrov, Tatjana

The expression of a gene is characterised by the upstream transcription factors and the biochemical reactions at the DNA processing them. Transient profile of gene expression then depends on the amount of involved transcription factors, and the scale of kinetic rates of regulatory reactions at the DNA. Due to the combinatorial explosion of the number of possible DNA configurations and uncertainty about the rates, a detailed mechanistic model is often difficult to analyse and even to write down. For this reason, modelling practice often abstracts away details such as the relative speed of rates of different reactions at the DNA, and how these reactions connect to one another. In this paper, we investigate how the transient gene expression depends on the topology and scale of the rates of reactions involving the DNA. We consider a generic example where a single protein is regulated through a number of arbitrarily connected DNA configurations, without feedback. In our first result, we analytically show that, if all switching rates are uniformly speeded up, then, as expected, the protein transient is faster and the noise is smaller. Our second result finds that, counter-intuitively, if all rates are fast but some more than others (two orders of magnitude vs. one order of magnitude), the opposite effect may emerge: time to equilibration is slower and protein noise increases. In particular, focusing on the case of a mechanism with four DNA states, we first illustrate the phenomenon numerically over concrete parameter instances. Then, we use singular perturbation analysis to systematically show that, in general, the fast chain with some rates even faster, reduces to a slow-switching chain. Our analysis has wide implications for quantitative modelling of gene regulation: it emphasises the importance of accounting for the network topology of regulation among DNA states, and the importance of accounting for different magnitudes of respective reaction rates. We conclude the paper by discussing the results in context of modelling general collective behaviour.

2022, Petrov, Tatjana, Tognazzi, Stefano

Spreading phenomena arise from simple local interaction among a large number of actors through different networks of interactions. Computational modelling and analysis of such phenomena is challenging due to the combinatorial explosion of possible network configurations. Traditional (single layer) networks are commonly used to encode the heterogeneous relationships among agents but are limited to a single type of interaction. Multiplex Multi-Layer networks (MLNs) have been introduced to allow the modeler to compactly and naturally describe multiple types of interactions and multiple simultaneous spreading phenomena. The downside is an increase in the complexity of the already challenging task of the analysis and simulation of such spreading processes. In this paper we explore the use of lumping techniques that preserve dynamics, previously applied to Continuous Time Markov Chains (CTMC) and single layer networks to multiple spreading processes on MLNs.

2021-10-06, Petrov, Tatjana, Tognazzi, Stefano

The concept of role equivalence has been applied in social network analysis for decades. Early definitions recognized two social actors as role equivalent, if they have identical relationships to the same other actors. Although this rather strong equivalence requirement has been relaxed in different ways, it is often challenging to detect interesting, non-trivial role equivalences, especially for social networks derived from empirical data. Multi-layer networks (MLNs) are increasingly gaining popularity for modelling collective adaptive systems, for example, engineered cyber-physical systems or animal collectives. Multiplex networks, a special case of MLNs, transparently and compactly describe such complex interactions (social, biological, transportation), where nodes can be connected by links of different types. In this work, we first propose a novel notion of exact and approximate role equivalence for multiplex MLNs. Then, we implement and experimentally evaluate the algorithm on a suite of real-world case studies. Results demonstrate that our notion of approximate role assignment not only obtains non-trivial partitions over nodes and layers as well, but it provides a fine-grained hierarchy of role equivalences, which is impossible to obtain by (combining) the existing role detection techniques. We demonstrate the latter by interpreting in detail the case study of Florence families, a classical benchmark from literature.

2020-09, Beica, Andreea, Feret, Jérôme, Petrov, Tatjana

Biochemical molecules interact through modification and binding reactions, giving raise to a combinatorial number of possible biochemical species. The time-dependent evolution of concentrations of the species is commonly described by a system of coupled ordinary differential equations (ODEs). However, the analysis of such high-dimensional, non-linear system of equations is often computationally expensive and even prohibitive in practice. The major challenge towards reducing such models is providing the guarantees as to how the solution of the reduced model relates to that of the original model, while avoiding to solve the original model.

In this paper, we have designed and tested an approximation method for ODE models of biochemical reaction systems, in which the guarantees are our major requirement. Borrowing from tropical analysis techniques, we look at the dominance relations among terms of each species' ODE. These dominance relations can be exploited to simplify the original model, by neglecting the dominated terms. As the dominant subsystems can change during the system's dynamics, depending on which species dominate the others, several possible modes exist. Thus, simpler models consisting of only the dominant subsystems can be assembled into hybrid, piecewise smooth models, which approximate the behavior of the initial system. By combining the detection of dominated terms with symbolic bounds propagation, we show how to approximate the original model by an assembly of simpler models, consisting in ordinary differential equations that provide time-dependent lower and upper bounds for the concentrations of the initial model's species.

The utility of our method is twofold. On the one hand, it provides a reduction heuristics that performs without any prior knowledge of the initial system's behavior (i.e., no simulation of the initial system is needed in order to reduce it). On the other hand, our method provides sound interval bounds for each species, and hence can serve to evaluate the faithfulness of tropicalization reduction heuristics for ODE models of biochemical reduction systems. The method is tested on several case studies.