What do molecules do when we are not looking? State sequence analysis for stochastic chemical systems

Abstract

Many biomolecular systems depend on orderly sequences of chemical transformations or reactions. Yet, the dynamics of single molecules or small-copy-number molecular systems are significantly stochastic. Here, we propose state sequence analysis--a new approach for predicting or visualizing the behaviour of stochastic molecular systems by computing maximum probability state sequences, based on initial conditions or boundary conditions. We demonstrate this approach by analysing the acquisition of drug-resistance mutations in the human immunodeficiency virus genome, which depends on rare events occurring on the time scale of years, and the stochastic opening and closing behaviour of a single sodium ion channel, which occurs on the time scale of milliseconds. In both cases, we find that our approach yields novel insights into the stochastic dynamical behaviour of these systems, including insights that are not correctly reproduced in standard time-discretization approaches to trajectory analysis.

Figure 1.

Alternative approaches to analysing probable pathwise behaviours of continuous-time Markov chains. (a) An example chain with three states. Dwell time parameters and transition probabilities are shown. (b) Ten randomly sampled trajectories of the system, starting from state 1, for a period of 6 s. Each trajectory comprises a specific sequence of states and dwell times. (c) The probabilities that the system transits different sequences of states as a function of time, after averaging over the possible transition times. The coloured bar at the top indicates the single most probable state sequence as a function of time, obtained using state sequence analysis. (d) A discrete-time Markov chain approximating the continuous-time chain. (e) Maximum probability trajectories for the time-discretized chain, obtained using standard dynamic programming techniques. Results differ from (c) because (e) represents single trajectories, whereas the sequences analysed in (c) represent integration over all possible trajectories having the same state sequence. (f) Stochastic simulation can be used to estimate state sequence probabilities, but there are uncertainties in the estimates, which sometimes results in incorrect identification of the maximally probable sequence.

Figure 2.

Analysis of within-patient mutational dynamics of HIV subject to Efavirenz combination therapy. (a) State transition diagram of the model, estimated based on time-series HIV genotype observations [34]. States are labelled with mutations they include. Details of the estimation method and all fitted model parameters can be found in §4.5. States participating in maximum probability state sequences are highlighted in yellow. (b) Probabilities of the most probable state sequences, starting from wild-type, as a function of time during the first 10 years of therapy. Coloured bars along the bottom show the time intervals during which different state sequences (indicated by colour) are maximally probable, according to state sequence analysis, stochastic simulation and time-discretized dynamic programming. (c) Depiction of the 100 most probable state sequences at different times. Each row corresponds to a different path, with the sequence of coloured rectangles depicting different states according to the legend at the right. The most probable path is on the top row, and the 100th most probable path is on the bottom row.

Figure 3.

Analysis of the dynamics of a single neural sodium channel. (a) Diagram of the model proposed in Vandenberg & Bezanilla [53], describing the gating dynamics of a single sodium channel in the squid giant axon. Transition rates a, b, c, d, f, g, i and j depend on the voltage of the clamp, as detailed in Vandenberg & Bezanilla [53] (see also §4.6). (b) Illustration of the relationship between current recordings of a single ion channel maintained at constant voltage via patch clamp (top), and an associated state trajectory simulated from the model in (a). Green and red coloured bars indicate open and closed periods, respectively. (c) Maximum probability sequences of closed states depending on the patch clamp voltage and the duration that the channel remains closed, obtained via state sequence analysis. Colour key to state sequences is below the panel. (d) Maximum probability sequences of closed states obtained by a time-discretization and dynamic programming. Colour key is below the panel.