Robustness and Games Against Nature in Molecular Programming

TL;DR

This paper introduces a game-theoretic method to assess the robustness of molecular programs, especially DNA-based ones, against perturbations in reaction timing, crucial for safety-critical nanoscale applications.

Contribution

It presents a novel game-theoretic approach to quantify molecular program robustness against clock speed perturbations, including random and malicious manipulations.

Findings

The method effectively measures robustness against random perturbations.

It provides a framework for assessing security against malicious manipulations.

Preliminary results demonstrate the approach's potential for safety-critical applications.

Abstract

Matter, especially DNA, is now programmed to carry out useful processes at the nanoscale. As these programs and processes become more complex and their envisioned safety-critical applications approach deployment, it is essential to develop methods for engineering trustworthiness into molecular programs. Some of this can be achieved by adapting existing software engineering methods, but molecular programming also presents new challenges that will require new methods. This paper presents a method for dealing with one such challenge, namely, the difficulty of ascertaining how robust a molecular program is to perturbations of the relative "clock speeds" of its various reactions. The method proposed here is game-theoretic. The robustness of a molecular program is quantified in terms of its ability to win (achieve its original objective) in games against other molecular programs that manipulate its relative clock speeds. This game-theoretic approach is general enough to quantify the security of a molecular program against malicious manipulations of its relative clock speeds. However, this preliminary report focuses on games against nature, games in which the molecular program's opponent perturbs clock speeds randomly (indifferently) according to the probabilities inherent in chemical kinetics.

Equations

$$\rho =(\mathbf{r},\mathbf{p},k)\in {\mathbb{N}}^S\times {\mathbb{N}}^S\times (0,\infty ),$$

$$X+C\stackrel{k}{\to }2Y+C$$

$$\mathrm{arity}(\rho )=\sum _{Y\in S}\mathbf{r}(\rho )(Y),$$

$$3Y+Z\to \text{RHS}.$$

$${\mathrm{rate}}_{\mathbf{x}}(\rho )$$

$$=ky(y-1)(y-2)z/V^3.$$

$$\mathrm{rate}(\mathbf{x},\mathbf{y})=\sum _{\rho (\mathbf{x})=\mathbf{y}}{\mathrm{rate}}_{\mathbf{x}}(\rho ).$$

$$U_i(\sigma ,\xi )={\int }_{{\mathbb{T}}_{\mathcal{G}}}{u}_i(\tau ){\mu }_{\sigma ,\xi }(\tau ),$$

$$U_1((N_1,N_2,\dots ,N_n),\xi )\ge \alpha U_1((N_1,(\varnothing ,\varnothing ),\dots ,(\varnothing ,\varnothing )),{\xi }^{\prime }),$$

$$R:\begin{array}[]{l}2X+Y\overset{1}{\to}3X\\ X+2Y\overset{1}{\to}3Y\end{array}$$

$$R^{\prime}:\begin{array}[]{l}2X+Y+A\overset{1}{\to}3X+A\\ X+2Y+B\overset{1}{\to}3Y+B.\end{array}$$

$$N_{k}:\begin{array}[]{l}A\overset{k}{\to}B\\ B\overset{k}{\to}A\end{array},$$

Full text

1 Introduction

Molecular programming is, at its simplest, computation with DNA. A programmed molecular system is a nanosystem that will execute the algorithmic behavior encoded into it. Examples of programmed DNA systems include neural net simulation, probabilistic switching circuits, nano-robotic walkers, and oscillators [21, 27, 7, 24]. That is, we are programming matter itself when we create a programmed molecular system.

Many of the intended uses of molecular programming are safety-critical, such as biosensors to detect pollutants, diagnostic devices to identify diseases, and nano-robotic walkers to perform targeted drug delivery [9]. Software engineering techniques including goal-oriented requirements modeling [26], risk analyses [14], and probabilistic model checking [15] have recently been extended to the emerging field of molecular programming, in order to aid the development of safe programmed molecular systems [16, 15, 10].

Assuring the robustness of a molecular program needs to occur before such a system is deployed. However, what robustness means for such a system is not well-defined. This, in turn, hinders efforts to determine how robust a particular system is.

The problem addressed by this paper is the difficulty of answering the question, “how robust is this molecular program?” The contribution of the paper is to propose a game-theoretic method by which we can quantitatively evaluate how robust a molecular program is to an opponent’s perturbing the relative clock speeds of its constituent processes (reactions). This sort of robustness is especially important, because the “rate constants” that govern the rates of chemical reactions are notoriously approximate and non-constant in actual laboratory experiments. We formulate this as a game against nature [17, 19], in which nature manipulates clock speeds at random, disinterested in the outcome. Although this approach is general enough to evaluate robustness in the face of a game against a malicious opponent, we focus here on the random perturbations inflicted by an indifferent nature. We illustrate our method on an important consensus algorithm, approximate majority [5]. We thus develop a games-against-nature formalism of robustness and then evaluate it on a small molecular program that computes the approximate majority.

2 Molecular Programs

Molecular programs are typically specified as chemical reaction networks (CRNs) [1, 6], which are roughly equivalent to stochastic Petri nets [8]. These CRNs can then be automatically compiled into DNA strand displacement systems that can be implemented in laboratory experiments [23, 4, 2].

Syntax

We now review the definition of CRNs. We fix a countably infinite set $\mathbf{S}$ whose elements are called species. We informally regard each species as an abstract type of molecule.

A reaction over a finite set $S\subseteq\mathbf{S}$ is an ordered triple

[TABLE]

where $\mathbf{r}\neq\mathbf{p}$ and $\mathbb{N}^{S}$ is the set of functions from $S$ into $\mathbb{N}$. Given such a reaction $\rho$, we call $\mathbf{r}(\rho)=\mathbf{r}$ the reactant vector of $\rho$, $\mathbf{p}(\rho)=\mathbf{p}$ the product vector of $\rho$, and $k(\rho)=k$ the rate constant of $\rho$. (Since $S$ is finite, it is natural to regard elements of $\mathbb{N}^{S}$ as vectors.) The species in the support $\operatorname{supp}(\mathbf{r})=\{X\in S\mid\mathbf{r}(X)>0\}$ are the reactants of $\rho$, and the species in $\operatorname{supp}(\mathbf{p})$ are the products of $\rho$.

We usually write reactions in a more chemistry-like notation. For example, if $S=\{X,Y,C\}$, then we write

[TABLE]

for the reaction $(\mathbf{r},\mathbf{p},k)$, where $\mathbf{r}(X)=1$, $\mathbf{r}(Y)=0$, $\mathbf{r}(C)=1$, $\mathbf{p}(X)=0$, $\mathbf{p}(Y)=2$, and $\mathbf{p}(C)=1$. A species $C$ satisfying $\mathbf{r}(\rho)(C)=\mathbf{p}(\rho)(C)>0$, as in this example, is called a catalyst of the reaction $\rho$. The net effect of a reaction $\rho$ is the vector $\Delta\rho=\mathbf{p}(\rho)-\mathbf{r}(\rho)\in\mathbb{Z}^{S}$. The arity of a reaction $\rho$ is

[TABLE]

i.e., its total number of reactants. A chemical reaction network (CRN) is an ordered pair $N=(S,R)$, where $S\subseteq\mathbf{S}$ is finite and $R$ is a finite set of reactions over $S$.

Semantics

In this paper we assign each CRN $N=(S,R)$ the operational meaning given by the stochastic mass action semantics (also called the stochastic mass action kinetics) introduced by Gillespie [13]. In this semantics a state of $N$ is a vector $\mathbf{x}\in\mathbb{N}^{S}$. For each $Y\in S$, the component $\mathbf{x}(Y)$ of $\mathbf{x}$ is the count of species $Y$ in the state $\mathbf{x}$. A reaction $\rho$ is applicable to state $\mathbf{x}$ if $\mathbf{r}(\rho)\leq\mathbf{x}$, i.e., all the reactants of $\rho$ are present in $\mathbf{x}$. If $\rho$ is applicable to $\mathbf{x}$, then the result of applying $\rho$ to $\mathbf{x}$ is the state $\rho(\mathbf{x})=\mathbf{x}+\Delta\rho\in\mathbb{N}^{S}$.

The (stochastic mass action) rate of a reaction $\rho$ in a state $\mathbf{x}\in\mathbb{N}^{S}$ and volume $V>0$ of solution, which we denote by $\operatorname{rate}_{\mathbf{x}}(\rho)$, was defined and justified by Gillespie [13]. Here we give a single example. Let $\rho$ be a reaction

[TABLE]

(The right-hand side RHS does not affect the rate of a reaction.) For brevity, write $y=\mathbf{x}(Y)$ and $z=\mathbf{x}(Z)$. If $\rho$ is applicable to $\mathbf{x}$ (i.e., if $y\geq 3$ and $z\geq 1$), then

[TABLE]

Under stochastic mass-action semantics, a CRN $N=(S,R)$ functions as a continuous-time Markov chain [22] with state space $\mathbb{N}^{S}$ and, for each $\mathbf{x},\mathbf{y}\in\mathbb{N}^{S}$, transition rate

[TABLE]

The CRN $N$ is initialized to some state or distribution over states. When it enters a state $\mathbf{x}$, it stays there for a random, real-valued sojourn time $t\in(0,\infty]$ before instantaneously executing some reaction $\rho$ and jumping to the state $\rho(\mathbf{x})$. A trajectory of $N$ is thus a sequence $\tau=((\mathbf{x}_{0},t_{0}),(\mathbf{x}_{1},t_{1}),\ldots)$ of ordered pairs $(\mathbf{x}_{i},t_{i})$, where $\mathbf{x}_{i}$ is a state of $N$ and $t_{i}$ is the associated sojourn time. The trajectory $\tau$ is finite if it reaches a state to which no reaction is applicable. Otherwise, the trajectory is infinite.

3 CRN Games

To quantify the robustness of a CRN, we consider how its performance might be affected by other CRNs that are present in the same solution. Clearly, this evaluation will depend both on how performance is defined and on what kinds of other CRNs are present. We begin by describing a general game-theoretic framework that allows for any scalar quantification of performance (by defining appropriate utility functions) and arbitrary constraints placed on the other CRNs (by restricting the other players’ strategy spaces). We then discuss the special case where interactions between CRNs are mediated only by catalysts.

An $n$-player CRN game with players $1,2,\ldots,n$ is a pair $\mathcal{G}=(\mathcal{N},\mathbf{u})$ with the following components.

$\mathcal{N}=\mathbf{N}_{1}\times\mathbf{N}_{2}\times\ldots\times\mathbf{N}_{n}$ is the strategy profile space. To play the game, each player $i$ selects a strategy: a CRN $N_{i}=(S_{i},R_{i})$ from its strategy space $\mathbf{N}_{i}$, which is a set of CRNs. For convenience we require that $(\emptyset,\emptyset)\in\mathbf{N}_{i}$. The $n$ players’ selected strategies collectively define a strategy profile $\sigma\in\mathcal{N}$ and a CRN $N_{\sigma}=(S_{\sigma},R_{\sigma})$, where $S_{\sigma}=\bigcup_{i}S_{i}$ and $R_{\sigma}=\bigcup_{i}R_{i}$. A state of the game comprises counts of all species in $S_{\sigma}$. As in a CRN, a trajectory in this game is a (finite or infinite) sequence of states paired with sojourn times. The space of all trajectories for the game $\mathcal{G}$ is $\mathbb{T}_{\mathcal{G}}$. 2. 2.

$\mathbf{u}=(u_{1},u_{2},\ldots,u_{n})$ is a profile of utility functions $u_{i}:\mathbb{T}_{\mathcal{G}}\to\mathbb{R}$, where $u_{i}(\tau)$ is the utility player $i$ gets from trajectory $\tau$. In our example below, this utility function is simply a binary indicator for “success,” meaning that the player gets utility 1 from any trajectory that performs a given task correctly and gets utility 0 from all other trajectories. A player representing “nature” is totally indifferent to the game’s outcome and hence gets utility 0 from all trajectories.

Let $\xi=(\xi_{1},\xi_{2},\ldots,\xi_{n})$ where $\xi_{i}\in\mathbb{N}^{S_{i}}$ is the random vector for $N_{i}$’s initial state, and let $\hat{\xi}_{i}$ be the embedding of $\xi_{i}$ in $\mathbb{N}^{S_{\sigma}}$. Then the initial state of $N_{\sigma}$ is the random vector $\sum_{i}\hat{\xi}_{i}$. Given this initial distribution, the theory of continuous-time Markov chains specifies a probability measure on the set $\mathbb{T}_{N_{\sigma}}$ of all trajectories of the CRN $N_{\sigma}$ [22], which also immediately yields a probability measure $\mu_{\sigma,\xi}$ on $\mathbb{T}_{\mathcal{G}}$. This allows us to define the function

[TABLE]

which is the expected value of player $i$’s utility when the strategy profile $\sigma$ is played.

Robustness

We measure the robustness of a CRN $N_{1}$ to a profile $(N_{2},N_{3},\ldots,N_{n})$ of other players’ CRNs by comparing player 1’s expected utility playing $N_{1}$ against those CRNs to its expected utility playing $N_{1}$ against trivial CRNs. Formally, for any $\alpha\in[0,1]$, a CRN $N_{1}\in\mathbf{N}_{1}$ is $\alpha$-robust to $(N_{2},N_{3},\ldots,N_{n})$ in game $\mathcal{G}$ under $\xi$ if

[TABLE]

where $\xi^{\prime}=(\xi_{1},\epsilon,\epsilon,\ldots,\epsilon)$ and $\epsilon$ is the trivial distribution that assigns probability 1 to the empty vector. This means that the participation of other players using these strategies can decrease player 1’s expected utility by at most a factor of $\alpha$.

Catalytic Games

The very general CRN games that we have defined allow essentially unrestricted interactions among the players’ CRNs. For many purposes, including those of this paper, it is more appropriate to restrict these interactions to those mediated by catalysts.

A catalytic game is one in which each player’s set of species can be written as $\mathbf{S}_{i}=\mathbf{A}_{i}\cup\mathbf{C}_{i}$ such that

$\mathbf{A}_{i}$ and $\mathbf{A}_{j}$ are disjoint for all $j\neq i$, and 2. 2.

for all $C\in\mathbf{C}_{i}$ and $\rho\in\mathbf{R}_{i}$, we have $\Delta\rho(C)=0$.

In such a game, each player can affect other players only by altering the counts of their catalysts, and hence only by altering the rates of their reactions.

4 Example

Approximate Majority

In this preliminary report, we use a game against nature to investigate the robustness of a simple chemical reaction network that computes approximate majority.

The task in approximate majority is to design a chemical reaction network $N$ with two designated species $X$ and $Y$ and the following objective. Let $x(t)$ and $y(t)$ be the counts of $X$ and $Y$, respectively, at time $t$. First, the total population $x(t)+y(t)$ should be constant as $t$ varies. Moreover, if $x(0)$ and $y(0)$ differ by a non-negligible amount, then whichever is larger should eventually “take over.” That is, if $x(0)\gg y(0)$, then we should with high probability have $x(t)=x(0)+y(0)$ (and hence $y(t)=0$) for all sufficiently large $t$. Similarly, if $y(0)\gg x(0)$, then we should with high probability have $y(t)=x(0)+y(0)$ for all sufficiently large $t$. If $x(0)$ is very close to $y(0)$, then we want one of these “takeovers” to occur with high probability, but it may be either species that takes over.

We investigate the robustness of the approximate majority CRN

[TABLE]

of Condon, Hajiaghayi, Kirkpatrick and Manuch [5]. In order to do this in a catalytic game, we replace $R$ with the catalyzed CRN

[TABLE]

The crucial thing to note here is that, if $x$, $y$, $a$, and $b$ are the counts of $X$, $Y$, $A$, and $B$ at some time, then the rates (“clock speeds”) of the reactions in $R$ at this time are $x(x-1)y$ and $xy(y-1)$, respectively, while the rates of the reactions in $R^{\prime}$ are $ax(x-1)y$ and $bxy(y-1)$, respectively. If $a$ and $b$ are equal and constant, then $R^{\prime}$ is merely a uniformly sped-up version of $R$. However, if $a$ and $b$ vary randomly, then the relative rates of the reactions in $R^{\prime}$ also vary randomly (i.e., the ratio of these rates varies randomly).

In order to test the robustness of $R$ to random perturbations of the relative rates of their reactions we thus play the CRN $R^{\prime}$ against a random “nature” that varies the initially equal counts $a$ and $b$ randomly. We model this behavior by the simple CRN

[TABLE]

which is “calibrated” by the rate constant $k\in(0,\infty)$.

Simulation

We assessed the robustness of the approximate majority CRN $R$ by comparing the performance of the catalyzed CRN $R^{\prime}$ in isolation with its performance in the presence of the CRN $N_{k}$ representing nature randomly perturbing rate constants. We frame this as a game where the utility of the approximate majority player is given by its success frequency.

We first created models in MATLAB using SimBiology software tools for the catalyzed approximate majority CRN $R^{\prime}$ and the nature CRN $N_{k}$ described above. With initial populations $a(0)=b(0)=100$ and combined initial population $x+y=10,000$, we varied the difference $x(0)-y(0)$ from 0 to 1,000 by intervals of 10. With these initial conditions, we ran $R^{\prime}$ in the presence of $N_{k}$ with $k=10^{9}$. We ran 10,000 trials for each set of initial conditions. Each trial converged within $10^{-8}$ time units, meaning that either $x(10^{-8})$ or $y(10^{-8})$ was 0. The design of the CRN $R^{\prime}$ guarantees that, once one population has taken over, no further reactions can occur.

As Fig. 1 shows, the random perturbations of $a$ and $b$ did reduce the success probability of the approximate majority algorithm for some values of $x(0)-y(0)$. For example, when $x(0)=5,120$ and $y(0)=4,880$, $R^{\prime}$ was $99\%$ successful in a vacuum but only $76\%$ successful in the presence of nature. However, even with the random perturbations, the success frequency in the presence of nature was always greater than $70\%$ of the success frequency in the absence of nature. This suggests that the CRN $R^{\prime}$ is at least $0.7$-robust to $N_{10^{9}}$ in this game under arbitrary distributions of initial states with $a=b=100$ and $x+y=10,000$.

5 Conclusion

Software engineering for molecular programming is a new research direction with open problems that can benefit from the attention of the software engineering research community. Many planned molecular systems will be deployed for use in vivo within a few years, and certification for safety-critical scenarios such as biosensors and drug delivery devices will require improved evidence of robustness. Molecular program developers similarly will be called upon to prevent system design vulnerabilities to malicious adversaries. Software engineering has an essential role to play in what scientists are already labeling as the century of life sciences [12].

The preliminary work described in this paper uses a game-theoretic approach to (1) formulate the robustness of a molecular program’s CRN model in terms of a game against nature and (2) provide a method to quantitatively evaluate its robustness. The example we present concerns random perturbations of the program’s clock speed by nature; however, the approach is general enough to also enable evaluation of security against an adversarial molecular program who maliciously perturbs the relative clock speeds. Future work on this will entail challenging issues involving strategic equilibria [18]. Our approach provides a foundation from which to pursue improved development and deployment of verifiably robust programmed molecular systems. More broadly, robustness in the presence of probabilistic behavior also is required for many non-molecular programmed systems [25, 11, 3, 20], and the advances described here may enhance our understanding of how to design in and verify robustness for other asynchronous systems operating in stochastic environments.

Acknowledgments

This research was supported in part by National Science Foundation grant 1545028. We thank Jim Lathrop for tool assistance.