inverse.tex

In this chapter
we would like to explore restricted versions of Kappa
for which it is possible to infer the energy function
from the rewriting rules and their associated rates.
Recall from the introduction
that in Kappa itself this problem is undecidable \citep{et1}.
One such restriction is when agents do not have sites
and thus cannot bind.
This is Petri nets.
We briefly present here the result obtained by \citet{et2}
and show the construction of the energy function
for \emph{simple} and \emph{symmetric} Petri nets
with \emph{mass action} semantics
(sisma Petri net for short).
\begin{definition}
A \emph{sisma Petri net} is a Petri net
on species $\species$ and reactions $\reactions$ for which:
\begin{enumerate}[label={(\roman*)}]
\item (simple) there are no two reactions that have
  the same stoichiometry (net change in species).
\item (symmetric) for each reaction $r \in \reactions$,
  there is a reaction $\inv{r} \in \reactions$
  that has the reverse direction,
  \ie the inputs of one are the outputs of the other.
\item (mass action) the jumping rate $q_{xy,r}$
  of going from a state $x$ to $y$
  by a reaction $r$ is proportional to
  the number of ocurrences of its left-hand side in $x$.
  In particular,
  given a rate constant $k(r)$ for reaction $r$,
  % the jumping rate of $r$ at state $x$ is
  we have
  \[ q_{xy,r} = k(r) \prod_{A \in \species}
     \frac{x(A)!}{(x(A)-\Delta_r(A))!} \]
  where $x(A)$ is the number of $A$s in $x$
  and $\Delta_r(A)$ is the net change of $A$ in reaction $r$.
\end{enumerate}
\end{definition}
Note that simple implies
no two distinct reactions can be applied
to a state $x$ to obtain state $y$.
A sisma Petri net has an energy function
\begin{equation}\label{eq:pn-energy}
  E(x) = \sum_{A \in \species} \cost(A) x(A) + \ln\bigl(x(A)!\bigr)
\end{equation}
for some function $\cost: \species \to \RR$ such that,
for all $r \in \reactions$,
\[ \sum_{A \in \species} \Delta_r(A) \cost(A) = \ln(k(r^\star)) - \ln(k(r)). \]
If there is no such function $\cost$,
the Petri net does not have detailed balance and
an energy function.
%
In the rest of the chapter
we introduce two other restrictions of Kappa
and show how to construct their energy function.


\section{Cooperative assembly systems}
\label{sec:cas}

The first restriction is when
rules can only create or destroy one edge at a time
and their rates can only depend on
how many bound sites the endpoints of the edge have.
% the type of the neighbours of the two endpoints of the edge
% but not the state of their sites
% (whether they are free or bound).
% Moreover, sites in an agent are regarded as \emph{indistinguishable}
% and thus rules can only count how many are bound.
Therefore sites are treated as \emph{indistinguishable}.
In addition, agents of the same type cannot bind.
\citet{cas} have proposed these restrictions
and a simple formalism incorporating them
to study the thermodynamics of polymer formation
when there are two types of monomers.
Here we extend their result to any number of monomer types.
% TODO: perhaps take the idea in the next sentence
% and integrate into the paragraph:
% Cooperativity means that the rate of binding (and unbinding)
% depends on the number of connections that have been established
% by the participating nodes.

In the case of two monomers,
rules are of the form
\begin{center}
  \begin{tikzpicture}
    \node[grphnode,anchor=east] (lhs) at (0,0) {
      \tikz[ingrphdiag]{
        \begin{scope}[shift={(0,0)}]
          \e{0,0}{90:1.1};
          \e{0,0}{150:1.1};
          \e{0,0}{270:1.1};
          % a
          \n[n1]{a}{0,0};
          \e{a}{.5,0};
          \site{x1}{a.north};
          \site{x2}{150:.25};
          \site{xn}{a.south};
          \site{xn+1}{a.east};
          % b1
          \n[n2]{b1}{90:1.1};
          \site{z}{b1.south};
          % b2
          \n[n2]{b2}{150:1.1};
          \site[shift={(150:1.1)}]{z}{-30:.25};
          % ellipsis
          \node[Black!60!White] at (200:.8) {\Large .};
          \node[Black!60!White] at (210:.8) {\Large .};
          \node[Black!60!White] at (220:.8) {\Large .};
          % bn
          \n[n2]{bn}{270:1.1};
          \site{z}{bn.north};
        \end{scope}
        \begin{scope}[shift={(1.2,0)}]
          \e{0,0}{90:1.1};
          \e{0,0}{30:1.1};
          \e{0,0}{270:1.1};
          % b
          \n[n2]{b}{0,0};
          \e{b}{-.5,0};
          \site{y1}{b.north};
          \site{y2}{30:7pt};
          \site{yn}{b.south};
          \site{yn+1}{b.west};
          % a1
          \n[n1]{a1}{90:1.1};
          \site{z}{a1.south};
          % a2
          \n[n1]{a2}{30:1.1};
          \site[shift={(30:1.1)}]{z}{210:.25};
          % ellipsis
          \node[Black!60!White] at (-20:.8) {\Large .};
          \node[Black!60!White] at (-30:.8) {\Large .};
          \node[Black!60!White] at (-40:.8) {\Large .};
          % an
          \n[n1]{an}{270:1.1};
          \site{z}{an.north};
        \end{scope}
      }};
    \path (lhs.east) +(.3,.09) edge[rule] +(1,.09)
      +(1.3,0) coordinate (r);
    \path (lhs.east) +(1,-.09) edge[rule] +(.3,-.09);
    \node[grphnode,anchor=west] (rhs) at (r) {
      \tikz[ingrphdiag]{
        \e{0,0}{1.1,0};
        \begin{scope}[shift={(0,0)}]
          \e{0,0}{90:1.1};
          \e{0,0}{150:1.1};
          \e{0,0}{270:1.1};
          % a
          \n[n1]{a}{0,0};
          \site{x1}{a.north};
          \site{x2}{150:.25};
          \site{xn}{a.south};
          \site{xn+1}{a.east};
          % b1
          \n[n2]{b1}{90:1.1};
          \site{z}{b1.south};
          % b2
          \n[n2]{b2}{150:1.1};
          \site[shift={(150:1.1)}]{z}{-30:.25};
          % ellipsis
          \node[Black!60!White] at (200:.8) {\Large .};
          \node[Black!60!White] at (210:.8) {\Large .};
          \node[Black!60!White] at (220:.8) {\Large .};
          % bn
          \n[n2]{bn}{270:1.1};
          \site{z}{bn.north};
        \end{scope}
        \begin{scope}[shift={(1.1,0)}]
          \e{0,0}{90:1.1};
          \e{0,0}{30:1.1};
          \e{0,0}{270:1.1};
          % b
          \n[n2]{b}{0,0};
          \site{y1}{b.north};
          \site{y2}{30:7pt};
          \site{yn}{b.south};
          \site{yn+1}{b.west};
          % a1
          \n[n1]{a1}{90:1.1};
          \site{z}{a1.south};
          % a2
          \n[n1]{a2}{30:1.1};
          \site[shift={(30:1.1)}]{z}{210:.25};
          % ellipsis
          \node[Black!60!White] at (-20:.8) {\Large .};
          \node[Black!60!White] at (-30:.8) {\Large .};
          \node[Black!60!White] at (-40:.8) {\Large .};
          % an
          \n[n1]{an}{270:1.1};
          \site{z}{an.north};
        \end{scope}
      }};
  \end{tikzpicture}
\end{center}
where the three grey dots on the sides of each graph
are an ellipsis to mean that monomers can be bound
to an arbitrary number of monomers of the other type
as long as each site is bound only once
and there is a finite number of sites per monomer
fixed by the monomer type.
% TODO: is this sentence clear?
Hence, the rule schema % formula?
represents a family of rules indexed by % $i,j$
the number of bound sites in the two monomers.
% TODO: add the following?
% Note that no nodes can be created or destroyed
% by this family of rules
% and so the total number of nodes is fixed.

We formalise the ideas in the first paragraph as follows.
Let $\types$ be the set of monomer types % $A,B,\dots$
and $\valence: \types \to \NN$ the map that assigns
% for each type the number of sites they have,
to each type the number of sites a monomer of that type has,
which we refer to as their valence.
A~monomer $u$ has type $\typeof(u) \in \types$
and degree $\degree_x(u) \in \NN$ in state $x$.
We simply write $\valence(u)$ for $\valence(\typeof(u))$.
The rate constant of a rule that binds
a monomer of type $\tp \in \types$ and degree $i$ with
a monomer of type $\tp' \in \types$ and degree $j$
is $\gamma^+_{\tp,i,\tp',j}$.
The rate constant of the reverse rule (unbinding)
is $\gamma^-_{\tp,i,\tp',j}$.
%
The jumping rate $q_{xy}$ from state $x$ to $y$
% is then determined by mass action semantics,
% \ie is linear in the number of ocurrences
% of the left-hand side of the rule that operates the transition.
is then linearly determined by the number of ocurrences
of the left-hand side of the rule and the rate constant
(mass action semantics).
We assume that any two agents can be bound only once.

The binding or unbinding of any two nodes $u,v$ in $x$
can only be carried out by one rule,
namely the one that operates on degrees $\degree_x(u),\degree_x(v)$.
The binding rule has rate constant
$\alpha(u,v) := \gamma^+_{\typeof(u),\degree_x(u),\typeof(v),\degree_x(v)}$
while the unbinding rate constant is
$\beta(u,v) := \gamma^-_{\typeof(u),\degree_x(u),\typeof(v),\degree_x(v)}$.
When binding we are free to choose
one site among the $\valence(u)-\degree_x(u)$ free sites of $u$
and one among the $\valence(v)-\degree_x(v)$ free sites of $v$
in order to apply the binding rule.
On the other hand, when unbinding we have only one choice,
namely removing the only edge between $u$ and $v$.
Hence, $q_{xy}$ is equal to $\alpha(u,v) \,
(\valence(u) - \degree_x(u)) \, (\valence(v) - \degree_x(v))$
% in the forward direction (binding)
when the binding rule is applied to $x$ to obtain $y$
and to $\beta(u,v)$
% in the backward direction (unbinding).
when unbinding.

The theorem below shows under which conditions
the type of systems presented in this section
have an energy function.

\begin{proposition}
  \label{prop:cas}
  Let $\types$ be a finite set of monomer types
  and $\gamma^-_{\tp,i,\tp',j},\gamma^+_{\tp,i,\tp',j}$
  families of real values indexed by types $\tp,\tp' \in \types$,
  $0 \leqslant i < \valence(\tp)$ and
  $0 \leqslant j < \valence(\tp')$ as above.
  Given a family $\Gamma_{\tp,i}$ of non-zero real values
  the following two statements are equivalent
  \begin{enumerate}[label={(\roman*)}]
  \item The \qmatrix $\qm$ as defined above by $q_{xy}$
    has detailed balance with respect to the \pmf $\pi$
    determined by the energy function
    \[ E(x) = \sum_{\tp \in \types} \sum_{0 < i \leqslant \valence(\tp)}
              \cost_\tp(i) \, x(\tp_i) \]
    where $x(\tp_i)$ is the number of nodes of type $\tp$
    with degree $i$ in $x$ and
    % TODO: which other variable name is good for the iterator a?
    \[ \cost_\tp(i) = \sum_{0 \leqslant j < i}
       \ln\, \frac{\Gamma_{\tp,j}}{\valence(\tp) - j} \]
  % \item The ratio $\gamma^-_{\tp,i,\tp',j} / \gamma^+_{\tp,i,\tp',j}$
  %   is equal to $\Gamma_{\tp,i} \, \Gamma_{\tp',j}$
  \item For all $\tp,\tp' \in \types$,
    $0 \leqslant i < \valence(\tp)$ and
    $0 \leqslant j < \valence(\tp')$ we have
    \begin{equation}
      \label{eq:cas-ratio}
      \frac{\gamma^-_{\tp,i,\tp',j}}{\gamma^+_{\tp,i,\tp',j}} =
      \Gamma_{\tp,i} \, \Gamma_{\tp',j}
    \end{equation}
  \end{enumerate}
\end{proposition}
\begin{proof}
  ($\Rightarrow$):
  Recall the detailed balance condition
  from \defn{detailed-balance}
  which says that for all states $x,y$
  \[ \pi_x \, q_{xy} = \pi_y \, q_{yx} \]
  By substituting $\pi_x$ and $\pi_y$ as in \eqn{energy} we obtain
  \begin{equation*}
    \expn\left[
      -\sum_{\tp \in \types} \sum_{0 < i \leqslant \valence(\tp)}
      \cost_\tp(i) \, x(\tp_i)
    \right] q_{xy} =
    \expn\left[
      -\sum_{\tp \in \types} \sum_{0 < i \leqslant \valence(\tp)}
      \cost_\tp(i) \, y(\tp_i)
    \right] q_{yx}
  \end{equation*}
  and by rearranging
  \[ \prod_{\tp \in \types} \prod_{0 < i \leqslant \valence(\tp)}
     e^{\,\cost_\tp(i) \, (y(\tp_i) - x(\tp_i))} =
     \frac{q_{yx}}{q_{xy}} \]

  When $y$ is obtained from $x$ by binding nodes $u,v$,
  the difference $y(\tp_i) - x(\tp_i)$ is equal to $0$
  for all pairs $\tp, i$ except
  i) when $\tp = \typeof(u), i = \degree_x(u)$
  or $\tp = \typeof(v)$, $i = \degree_x(v)$,
  then $y(\tp_i) - x(\tp_i) = -1$; and
  ii) when $\tp = \typeof(u), i = \degree_y(u) = \degree_x(u) + 1$
  or $\tp = \typeof(v), i = \degree_y(v) = \degree_x(v) + 1$,
  then $y(\tp_i) - x(\tp_i) = 1$.
  Let $\tp_u = \typeof(u)$, $\tp_v = \typeof(v)$,
  $d_u = \degree_x(u)$ and $d_v = \degree_x(v)$.
  It follows that the last equation can be rewritten as
  \[ \expn\left[
       \cost_{\tp_u}(d_u+1) +
       \cost_{\tp_v}(d_v+1) -
       \cost_{\tp_u}(d_u) -
       \cost_{\tp_v}(d_v) \right] =
     \frac{q_{yx}}{q_{xy}} \]
  % By expanding ... % Philipp didn't like expanding
  By substituting $\cost$ we get % and $q$ we have
  % \[ \frac{
  %      \Gamma_{\tp_u,d_u} \, \Gamma_{\tp_v,d_v}}{
  %      (\valence(u) - d_u) \, (\valence(v) - d_v)} =
  \[ \frac{
     \prod_{0 \leqslant i < d_u+1}%\limits
     \frac{\Gamma_{\tp_u,i}}{\valence(u) - i} \,
     \prod_{0 \leqslant i < d_v+1}%\limits
     \frac{\Gamma_{\tp_v,i}}{\valence(v) - i}}{
     \prod_{0 \leqslant i < d_u}%\limits
     \frac{\Gamma_{\tp_u,i}}{\valence(u) - i} \,
     \prod_{0 \leqslant i < d_v}%\limits
     \frac{\Gamma_{\tp_v,i}}{\valence(v) - i}} =
     \frac{q_{yx}}{q_{xy}} \]
     % \frac{
     %   % \gamma^-_{\tp_u,d_u,\tp_v,d_v}}{
     %   % \gamma^+_{\tp_u,d_u,\tp_v,d_v} \,
     %   \beta(u,v)}{
     %   \alpha(u,v) \, (\valence(u) - d_u) \, (\valence(v) - d_v)} \]
  Products on the left cancel out and yield,
  after substituting $q$ on the right,
  % which by simplification and substitution of $q$ yields
  \[ \frac{
       \Gamma_{\tp_u,d_u} \, \Gamma_{\tp_v,d_v}}{
       (\valence(u) - d_u) \, (\valence(v) - d_v)} =
     \frac{
       \beta(u,v)}{
       \alpha(u,v) \, (\valence(u) - d_u) \, (\valence(v) - d_v)} \]
  which then simplifies to
  \[ \Gamma_{\tp_u,d_u} \, \Gamma_{\tp_v,d_v} =
    \frac{\gamma^-_{\tp_u,d_u,\tp_v,d_v}}{\gamma^+_{\tp_u,d_u,\tp_v,d_v}} \]
  This equality holds in general for nodes of any degree and type.

  ($\Leftarrow$):
  We prove that, whenever (ii) holds,
  $\pi$ verifies the detailed balance condition.
  For all $x,y$ such that $q_{xy} = 0$
  the equality $\pi_x\,q_{xy} = \pi_y\,q_{yx}$ holds
  as rules are reversible and (ii) dictates that
  a rate constant is zero if the reverse rate constant is.
  % NOTE: q_{xy} = 0 might be because there is no transition
  % between x and y or because the rate constant is zero.
  When $q_{xy} > 0$ then $y$ can be obtained from $x$
  by binding or unbinding some nodes $u,v$.
  By substituting $\tp$ for $\typeof(u)$, $\tp'$ for $\typeof(v)$,
  $i$ for $\degree_x(u)$ and $j$ for $\degree_x(v)$
  in \eqn{cas-ratio} we obtain the last equation
  in the first part of the proof.
  We can replay the transformations backwards
  to obtain $\pi_x\,q_{xy} = \pi_y\,q_{yx}$
  when $y$ is obtained by binding.
  The case of unbinding follows a similar argument.
\end{proof}


\section{Flipping and binding} % co-ANC
\label{sec:fb}

Now we look at systems whose nodes have sites
that possess an internal state.
This internal state is used to decide when to bind other nodes
or change the internal state of other sites.
For simplicity, internal states can
take one of only two possible values. %, 0 and 1.
Unlike cooperative assembly systems,
% where sites were undistinguishable from each other,
here sites are \emph{distinguishable}
% and thus have a unique name that distinguishes it from the rest
as in Kappa.
Hence, we extend contact maps $g$ as defined in~\sct{kappa}
with a map $\sitestate_g$
% that assign an internal state $\sitestate_g(i)$ in $\set{0,1}$
% to sites $i$ in $\sites_{\anon{g}}$
that assigns an internal state vector
$\sitestate_g(u)$ in $\set{0,1}^{\sitemap_C^{-1}(g_\agents(u))}$
to agents $u$ in $\agents_{\anon{g}}$.
The internal state vector is indexed by
the sites of the agent type $g_\agents(u)$ of $u$.
We use these extended contact maps as graphs in this section.
% TODO: add the order of sites wherever it's used.
% In addition, sites are \emph{ordered}
% by a total order $<_a$ with a an agent type
% (node in the contact graph).

A site's internal state can be changed by rules
we call \emph{flips}.
The rate at which we flip a site may depend on
the type of the site and the node it belongs to,
the internal state of sites on the same node
and the type of the neighbours.
Note that it cannot depend on
the internal states of the neighbours' sites
or the nodes they are bound to.
Graphically, flips are of the form
\begin{center}
  \begin{tikzpicture}[thick]
    \node[grphnode,anchor=east] (lhs) at (0,0) {
      \tikz[ingrphdiag]{
        \nn[n4]{a}{0,0}{$u$};
        \n[n4,dotted]{b1}{10:1.1};
        \n[n4,dotted]{b2}{170:1.1};
        \n[n4,dotted]{b3}{-90:1.1};
        \e[dotted]{a}{b1};
        \e[dotted]{a}{b2};
        \e[dotted]{a}{b3};
        \site[n]{s1}{a.south};
        \site[n6]{s2}{170:.33};
        \site[n6]{s3}{10:.33};
        \site[n4,dotted,shift={(10:1.1)}]{s4}{190:.25};
        \site[n4,dotted,shift={(170:1.1)}]{s5}{-10:.25};
        \site[n4,dotted]{s6}{b3.north};
        \node at (-65:.48) {\scriptsize x};
        \node at (110:.7) {\Large .};
        \node at (90:.7) {\Large .};
        \node at (70:.7) {\Large .};
      }};
    \path (lhs.east) +(.3,.09) edge[rule] +(1,.09)
      +(1.3,0) coordinate (r);
    \path (lhs.east) +(1,-.09) edge[rule] +(.3,-.09);
    \node[grphnode,anchor=west] (rhs) at (r) {
      \tikz[ingrphdiag]{
        \nn[n4]{a}{0,0}{$u$};
        \n[n4,dotted]{b1}{10:1.1};
        \n[n4,dotted]{b2}{170:1.1};
        \n[n4,dotted]{b3}{-90:1.1};
        \e[dotted]{a}{b1};
        \e[dotted]{a}{b2};
        \e[dotted]{a}{b3};
        \site[n5]{s1}{a.south};
        \site[n6]{s2}{170:.33};
        \site[n6]{s3}{10:.33};
        \site[n4,dotted,shift={(10:1.1)}]{s4}{190:.25};
        \site[n4,dotted,shift={(170:1.1)}]{s5}{-10:.25};
        \site[n4,dotted]{s6}{b3.north};
        \node at (-65:.48) {\scriptsize x};
        \node at (110:.7) {\Large .};
        \node at (90:.7) {\Large .};
        \node at (70:.7) {\Large .};
      }};
  \end{tikzpicture}
\end{center}
where the dotted lines % in the rule diagram
denote an optional node or edge
and site $x$ changes state from white to black.
As usual, we write $r_L$ for the left-hand side contact map of rule $r$
and $r_R$ for that of the right-hand side,
both contact maps over some fixed contact graph $C$.

In a flip we have a complete view
over the internal state of sites in $u$
and those of the neighbours,
which we characterise as vectors indexed by site types
in $I := \sitemap_C^{-1}(r_{L,\agents}(u))$.
The rate constants of the forward and backward
flip rules are then parametrised by the agent type
$a := r_{L,\agents}(u) = r_{R,\agents}(u)$ of $u$ in $C$,
the site type $i := r_{L,\sites}(x) = r_{R,\sites}(x)$ of $x$,
the internal state vector $\sitestates \in \set{0,1}^I$ of $u$,
% TODO: type of the internal state vector of the neighbour?
and the binding state vector $\neighbours \in
((\agents_C \times \sites_C) \union \set{\star})^I$
% ((\agents_C \times \sites_C \times \set{0,1}^I) \union \set{\star})^I
where $\star$ is used to denote a free site.
% When all sites are free we use $\emptyvec$.
We write $\lambda^+_{a,i,\sitestates,\neighbours}$
for the rate constant of the forward rule,
$\lambda^-_{a,i,\sitestates,\neighbours}$
for that of the backward rule,
and $\Lambda_{a,i,\sitestates,\neighbours} =
\lambda^-_{a,i,\sitestates,\neighbours}/
\lambda^+_{a,i,\sitestates,\neighbours}$ for their ratio.
% between the backward and forward rate constants. % of flips.
% TODO: do we (have to) use them?
% Also, later we use $\sitestatesof_x(u)$ and $\neighboursof_x(u)$
% for the internal state vector and neighbour vector % of neighbours' types
% of node $u$ in $x$.

The second type of rules we allow are \emph{binds}.
They are of the form
\begin{center}
  \begin{tikzpicture}[thick]
    \node[grphnode,anchor=east] (lhs) at (0,0) {
      \tikz[ingrphdiag]{
        \begin{scope}[shift={(0,0)}]
          \nn[n4]{a}{0,0}{$u$};
          \node[Black!60!White] at (155:.48) {\Large .};
          \node[Black!60!White] at (180:.48) {\Large .};
          \node[Black!60!White] at (205:.48) {\Large .};
          \node at (26:.48) {\scriptsize x};
          \e{a}{.55,0};
          \site[n6]{s1}{a.east};
          \site[n6]{s2}{110:.33};
          \site[n6]{s3}{-110:.33};
        \end{scope}
        \begin{scope}[shift={(1.4,0)}]
          \nn[n4]{b}{0,0}{$v$};
          \node[Black!60!White] at (25:.48) {\Large .};
          \node[Black!60!White] at (0:.48) {\Large .};
          \node[Black!60!White] at (-25:.48) {\Large .};
          \node at (206:.51) {\scriptsize y};
          \e{b}{-.55,0};
          \site[n6]{s1}{b.west};
          \site[n6]{s2}{70:.33};
          \site[n6]{s3}{-70:.33};
        \end{scope}
      }};
    \path (lhs.east) +(.3,.09) edge[rule] +(1,.09)
      +(1.3,0) coordinate (r);
    \path (lhs.east) +(1,-.09) edge[rule] +(.3,-.09);
    \node[grphnode,anchor=west] (rhs) at (r) {
      \tikz[ingrphdiag]{
        \begin{scope}[shift={(0,0)}]
          \nn[n4]{a}{0,0}{$u$};
          \node[Black!60!White] at (155:.48) {\Large .};
          \node[Black!60!White] at (180:.48) {\Large .};
          \node[Black!60!White] at (205:.48) {\Large .};
          \node at (26:.48) {\scriptsize x};
          \site[n6]{s2}{110:.33};
          \site[n6]{s3}{-110:.33};
        \end{scope}
        \begin{scope}[shift={(1.2,0)}]
          \nn[n4]{b}{0,0}{$v$};
          \node[Black!60!White] at (25:.48) {\Large .};
          \node[Black!60!White] at (0:.48) {\Large .};
          \node[Black!60!White] at (-25:.48) {\Large .};
          \node at (206:.51) {\scriptsize y};
          \site[n6]{s2}{70:.33};
          \site[n6]{s3}{-70:.33};
        \end{scope}
        \e{a}{b};
        \site[n6]{s1}{a.east};
        \site[n6]{s1}{b.west};
      }};
  \end{tikzpicture}
\end{center}
where the three grey dots are an ellipsis
meaning that nodes $u$ and $v$ must declare
an internal state to all sites they have.
The rate constant of the binding and unbinding rules
depends on the internal state of the sites
on the participating nodes $\sitestates$ and $\sitestates'$
as well as the type of the nodes $a,b$ and the bound sites~$i,j$.
We write $\Gamma_{a,i,\sitestates,b,j,\sitestates'} =
\gamma^-_{a,i,\sitestates,b,j,\sitestates'}/
\gamma^+_{a,i,\sitestates,b,j,\sitestates'}$ for
the ratio between the backward and forward rate constants of binds.

We refer to systems composed by flips and binds
as \emph{flip-bind systems} or FB-systems for short.
We will show how the energy function of an FB-system looks like.
But first, we prove two lemmas that show us how detailed balance
fixes the value of some ratios of rate constants,
reducing so the number of free parameters in the system.
To simplify notation in the following lemmas,
let $\sitestates+i$ be the vector $\sitestates$
with site $i$ flipped.

\begin{lemma}
  \label{lemma:flipflip}
  % Let $x$ be a state of an FB-system with detailed balance.
  % For all nodes $u$ in $x$ and sites $i,j$ in $u$ it is true that
  Let $C$ be a contact graph.
  For all agent types $a \in \agents_C$,
  let $I = \sitemap_C^{-1}(a)$,
  and for all
  site types $i \in I$, % \sites_C such that $\sitemap_C(i) = a$,
  internal state vectors
  $\sitestates \in \set{0,1}^I$ and
  binding state vectors
  $\neighbours \in ((\agents_C\times\sites_C)\union\set{\star})^I$,
  an FB-system with detailed balance verifies
  \begin{equation}
    \label{eq:flipflip}
    \Lambda_{a,i,\sitestates,\neighbours} \,
    \Lambda_{a,j,\sitestates+i,\neighbours} =
    \Lambda_{a,j,\sitestates,\neighbours} \,
    \Lambda_{a,i,\sitestates+j,\neighbours}
  \end{equation}
  % where $a$ is the type of $u$,
  % $\sitestates$ its internal state vector
  % and $\neighbours$ its vector of neighbours.
\end{lemma}
\begin{proof}
  Pick a state $x$ and a node $u$ in $x$.
  We can find a square in the transition graph
  with 4 flips starting from $x$:
  flip site $i$ first then $j$ and flip $j$ first then $i$.
  \begin{center}
    \begin{tikzpicture}
      \node (x) at (0,0) {$x$};
      \node (xi) at (4,0) {$x_i$};
      \node (xj) at (0,-3) {$x_j$};
      \node (xij) at (4,-3) {$x_{ij}$};
      % x -- xi
      \path (.5,.09) edge[rule]
        node[above] {$\lambda^+_{a,i,\sitestates,\neighbours}$}
        (3.5,.09);
      \path (3.5,-.09) edge[rule]
        node[below] {$\lambda^-_{a,i,\sitestates,\neighbours}$}
        (.5,-.09);
      % x -- xj
      \path (.09,-.5) edge[rule]
        node[right] {$\lambda^+_{a,j,\sitestates,\neighbours}$}
        (.09,-2.5);
      \path (-.09,-2.5) edge[rule]
        node[left] {$\lambda^-_{a,j,\sitestates,\neighbours}$}
        (-.09,-.5);
      % xi -- xij
      \path (4.09,-.5) edge[rule]
        node[right] {$\lambda^+_{a,j,\sitestates+i,\neighbours}$}
        (4.09,-2.5);
      \path (3.91,-2.5) edge[rule]
        node[left] {$\lambda^-_{a,j,\sitestates+i,\neighbours}$}
        (3.91,-.5);
      % x -- xi
      \path (.5,-2.91) edge[rule]
        node[above] {$\lambda^+_{a,i,\sitestates+j,\neighbours}$}
        (3.5,-2.91);
      \path (3.5,-3.09) edge[rule]
        node[below] {$\lambda^-_{a,i,\sitestates+j,\neighbours}$}
        (.5,-3.09);
    \end{tikzpicture}
  \end{center}
  By detailed balance we have that
  the product of rates along a cycle must be equal to 1.
  Hence, starting from $x$ and going through the cycle
  in one direction and the reverse we get
  \begin{equation*}
    \lambda^+_{a,i,\sitestates,\neighbours} \,
    \lambda^+_{a,j,\sitestates+i,\neighbours} \,
    \lambda^-_{a,i,\sitestates+j,\neighbours} \,
    \lambda^-_{a,j,\sitestates,\neighbours} =
    \lambda^+_{a,j,\sitestates,\neighbours} \,
    \lambda^+_{a,i,\sitestates+j,\neighbours} \,
    \lambda^-_{a,j,\sitestates+i,\neighbours} \,
    \lambda^-_{a,i,\sitestates,\neighbours}
  \end{equation*}
  By rearranging we obtain \eqn{flipflip}.
\end{proof}

When we consider all flips for a node with $n$ sites,
we find an $n$-hypercube in the transition graph.
This hypercube has $2^{n-1} n$ edges
where an edge corresponds to a pair of forward and backward flips.
It follows that there are the same number of $\Lambda$ ratios.
In addition, each face of the hypercube generates an equation
by \lem{flipflip} % detailed balance
and there are $2^{n-3} (n-1) n$ faces in an $n$-hypercube.
This is a severe constraint on the number of parameters
that can be freely set in an FB-system with detailed balance.
% Note also that the number of faces of an $n$-hypercube grows faster
% than the number of edges and thus when $n \geqslant 5$
% there are no free parameters.
% What does this mean practically?
% What is the concrete value that $\Lambda$s take?

A further constrain to the values of
the rate constant ratios of flips and bindings
can be obtained and its proved in the following lemma.
For this lemma we will use a fixed total order $<_a$
on the sites of an agent type $a$ in the contact graph $C$.

\begin{lemma}
  \label{lemma:flipbind}
  % Let $x$ be a state of an FB-system with detailed balance.
  % For all nodes $u$ in $x$ and sites $i$ in $u$ it is true that
  Let $C$ be a contact graph.
  For all agent types $a \in \agents_C$,
  let $I = \sitemap_C^{-1}(a)$,
  and for all
  site types $i \in I$, % \sites_C such that $\sitemap_C(i) = a$,
  binding state vectors
  $\neighbours \in ((\agents_C\times\sites_C)\union\set{\star})^I$
  and internal state vectors
  $\sitestates \in \set{0,1}^I$,
  % $\sitestates_{jb} \in \set{0,1}^{\sitemap_C^{-1}(b)}$
  % with $j \in I$, $b \in \agents_C$,
  $\sitestates_{j} \in \set{0,1}^{\sitemap_C^{-1}(b)}$
  with $j \in I$ and $b$ the agent type in $\neighbours(j)$,
  an FB-system with detailed balance verifies
  \begin{equation}
    \label{eq:flipbind}
    \frac{\Lambda_{a,i,\sitestates,\neighbours}}{
          \Lambda_{a,i,\sitestates,\emptyvec}} =
    \prod_{\substack{i' \in I\\\neighbours(i') \neq \star\\\tuple{b,j} := \neighbours(i')}}
    \frac{\Gamma_{a,i',\sitestates,b,j,\sitestates_{i'}}}{
          \Gamma_{a,i',\sitestates+i,b,j,\sitestates_{i'}}}
  \end{equation}
  % where $a$ is the type of $u$,
  % $\sitestates$ its internal state vector,
  % $\neighbours$ its vector of neighbours
  % and $\emptyvec$ a vector of free sites,
  where $\emptyvec$ a vector of free sites,
  \ie $\emptyvec(i) = \star$ for all $i \in I$.
\end{lemma}
\begin{proof}
  We construct a series of squares
  with 2 flips and 2 bindings each.
  Pick a state $x$ and a node $u$ in $x$.
  Strip $u$ of all its neighbours and
  call that state $x_0$.
  Choose a site $i$ of $u$.
  The first square starts from $x_0$, then
  i) flips site $i$ and
  ii) binds the smallest site $i_1$
  according to the order $<_{x_\agents(u)}$
  to site $j$ of an agent of type $b$
  that has internal vector state $\sitestates_{i_1}$,
  where $(j,b) = \neighbours(i_1)$.
  % or does (ii) first then (i).
  After performing (i) then (ii) or (ii) then (i)
  % the flip and bind in any of the two possible orders
  we reach state $x_1$.\\[.23cm]
  % \begin{center}
  % \vspace{.2cm}
  % \hspace*{-.3cm}
  \resizebox{\linewidth}{!}{%
    \begin{tikzpicture}
      %
      % first square
      %
      \node[outer sep=.2cm] (x0) at (0,0) {$x_0$};
      \node[outer sep=.2cm] (x0f) at (0,-3) {$\widehat{x_0}$};
      \node[outer sep=.2cm] (x1f) at (4,0) {$\widehat{x_1}$};
      \node[outer sep=.2cm] (x1) at (4,-3) {$x_1$};
      % x0 -- x1f
      \draw[rule,transform canvas={yshift=.09cm}] (x0) --
        node[above] {$\gamma^+_{a,i_1,\sitestates,b,j,\sitestates_{i_1}}$}
        (x1f);
      \draw[rule,transform canvas={yshift=-.09cm}] (x1f) --
        node[below] {$\gamma^-_{a,i_1,\sitestates,b,j,\sitestates_{i_1}}$}
        (x0);
      % x0 -- x0f
      \draw[rule,transform canvas={xshift=.09cm}] (x0) --
        node[right] {$\lambda^+_{a,i,\sitestates,\emptyvec}$}
        (x0f);
      \draw[rule,transform canvas={xshift=-.09cm}] (x0f) --
        node[left] {$\lambda^-_{a,i,\sitestates,\emptyvec}$}
        (x0);
      % x1f -- x1
      \draw[rule,transform canvas={xshift=.09cm}] (x1f) --
        node[right] {$\lambda^+_{a,i,\sitestates,\neighbours_1}$}
        (x1);
      \draw[rule,transform canvas={xshift=-.09cm}] (x1) --
        node[left] {$\lambda^-_{a,i,\sitestates,\neighbours_1}$}
        (x1f);
      % x0f -- x1
      \draw[rule,transform canvas={yshift=.09cm}] (x0f) --
        node[above] {$\gamma^+_{a,i_1,\sitestates+i,b,j,\sitestates_{i_1}}$}
        (x1);
      \draw[rule,transform canvas={yshift=-.09cm}] (x1) --
        node[below] {$\gamma^-_{a,i_1,\sitestates+i,b,j,\sitestates_{i_1}}$}
        (x0f);
      %
      % middle square
      %
      \node[outer sep=.2cm] (xn) at (7.5,0) {$x_n$};
      \node[outer sep=.2cm] (xnf) at (7.5,-3) {$\widehat{x_n}$};
      % x1f -- xn
      \draw[rule,transform canvas={yshift=.09cm}] (x1f) -- (xn);
      \draw[rule,transform canvas={yshift=-.09cm}] (xn) --
        node[below] {$\ldots$} (x1f);
      % xn -- xnf
      \draw[rule,transform canvas={xshift=.09cm}] (xn) --
        node[right] {$\lambda^+_{a,i,\sitestates,\neighbours_n}$}
        (xnf);
      \draw[rule,transform canvas={xshift=-.09cm}] (xnf) --
        node[left] {$\lambda^-_{a,i,\sitestates,\neighbours_n}$}
        (xn);
      % x1 -- xnf
      \draw[rule,transform canvas={yshift=.09cm}] (x1) --
        node[above] {$\ldots$} (xnf);
      \draw[rule,transform canvas={yshift=-.09cm}] (xnf) -- (x1);
      %
      % nth square
      %
      \node (xn+1) at (11.5,-3) {$x_{n+1}$};
      \node (xn+1f) at (11.5,0) {$\widehat{x_n}_{+1}$};
      % xn -- xn+1f
      \draw[rule,transform canvas={yshift=.09cm}] (xn) --
        node[above] {$\gamma^+_{a,i_n,\sitestates,b',j',\sitestates_{i_n}}$}
        (xn+1f);
      \draw[rule,transform canvas={yshift=-.09cm}] (xn+1f) --
        node[below] {$\gamma^-_{a,i_n,\sitestates,b',j',\sitestates_{i_n}}$}
        (xn);
      % xn+1f -- xn+1
      \draw[rule,transform canvas={xshift=.09cm}] (xn+1f) --
        node[right] {$\lambda^+_{a,i,\sitestates,\neighbours_{n+1}}$}
        (xn+1);
      \draw[rule,transform canvas={xshift=-.09cm}] (xn+1) --
        node[left] {$\lambda^-_{a,i,\sitestates,\neighbours_{n+1}}$}
        (xn+1f);
      % x0f -- x1
      \draw[rule,transform canvas={yshift=.09cm}] (xnf) --
        node[above] {$\gamma^+_{a,i_n,\sitestates+i,b',j',\sitestates_{i_n}}$}
        (xn+1);
      \draw[rule,transform canvas={yshift=-.09cm}] (xn+1) --
        node[below] {$\gamma^-_{a,i_n,\sitestates+i,b',j',\sitestates_{i_n}}$}
        (xnf);
    \end{tikzpicture}}
  % \end{center}
  % with $\neighbours_1$ the neighbour vector of $u$ in $x_1$.
  with $\neighbours_1(i) = \star$ for all $i$ except $i_1$,
  where $\neighbours_1(i_1) = \neighbours(i_1)$.
  In general, $\neighbours_n(i)$ is equal to $\neighbours(i)$
  if $i \leq_{x_\agents(u)} i_n$ and $\star$ otherwise.

  By detailed balance we obtain the following relation
  for the first square
  \begin{equation}
    \label{eq:fb1}
    \Lambda_{a,i,\sitestates,\emptyvec} \,
    \Gamma_{a,i_1,\sitestates+i,b,j,\sitestates_{i_1}} =
    \Gamma_{a,i_1,\sitestates,b,j,\sitestates_{i_1}} \,
    \Lambda_{a,i,\sitestates,\neighbours_1}
  \end{equation}

  We construct the $n$th square
  starting from $x_n$ by flipping site $i$
  and binding site $i_n$ to site $j'$ in an agent of type $b'$
  as indicated by $\neighbours(i_n)$.
  This neighbour has an internal state vector $\sitestates_{i_n}$.
  Again, by detailed balance we get
  \begin{equation}
    \label{eq:fbn}
    \Lambda_{a,i,\sitestates,\neighbours_n} \,
    \Gamma_{a,i,\sitestates+i,b',j',\sitestates_{i_n}} =
    \Gamma_{a,i,\sitestates,b',j',\sitestates_{i_n}} \,
    \Lambda_{a,i,\sitestates,\neighbours_{n+1}}
  \end{equation}
  % with $\neighbours_n$ the neighbour vector of $u$ in $x_n$.

  \eqn{fb1} can be rewritten as
  \begin{equation*}
    \frac{
      \Lambda_{a,i,\sitestates,\neighbours_1}}{
      \Lambda_{a,i,\sitestates,\emptyvec}} =
    \frac{
      \Gamma_{a,i,\sitestates+i,
              b,j,\sitestates_{i_1}}}{
      \Gamma_{a,i,\sitestates,
              b,j,\sitestates_{i_1}}}
  \end{equation*}
  Then we substitute $\Lambda_{a,i,\sitestates,\neighbours_1}$
  according to \eqn{fbn} and obtain
  \begin{equation*}
    \frac{
      \Lambda_{a,i,\sitestates,\neighbours_2}}{
      \Lambda_{a,i,\sitestates,\emptyvec}} =
    \frac{
      \Gamma_{a,i,\sitestates+i,
              b,j,\sitestates_{i_1}}}{
      \Gamma_{a,i,\sitestates,
              b,j,\sitestates_{i_1}}} \,
    \frac{
      \Gamma_{a,i,\sitestates+i,
              c,k,\sitestates_{i_2}}}{
      \Gamma_{a,i,\sitestates,
              c,k,\sitestates_{i_2}}}
  \end{equation*}
  We repeat until we recover \eqn{flipbind}.
\end{proof}

Now we proceed to compute the energy function.

\begin{proposition}
  \label{prop:fb-energy}
  Given an FB-system with detailed balance,
  its energy function is
  \begin{equation}
    \label{eq:fb-energy}
    % E(x) = \sum_{u \in \agents_{\anon{x}}
    % \sum_{s \in \sitemap_{\anon{x}}^{-1}(u)}
    E(x) = \sum_{s \in \sites_{\anon{x}}}
    \ln\,\Lambda_{a,i,\sitestates_u(i),\emptyvec} +
    \sum_{(s,t) \in \edges_{\anon{x}}}
    \ln\,\Gamma_{a,i,\sitestate_x(u),b,j,\sitestate_x(v)}
  \end{equation}
  where $u = \sitemap_{\anon{x}}(s)$, $v = \sitemap_{\anon{x}}(t)$,
  $a = x_\agents(u)$, $b = x_\agents(v)$,
  $i = x_\sites(s)$, $j = x_\sites(t)$,
  and $\sitestates_w: \sites_C \to \set{0,1}^I$
  a family of functions indexed by $w \in \agents_{\anon{x}}$,
  with $I = \sitemap_C^{-1}(x_\agents(w))$,
  defined as
  \begin{equation*}
    \sitestates_u(i)(j) = \begin{cases}
      \sitestate_x(u)(j) & \text{if $j <_{x_\agents(u)} i$} \\
      0 & \text{otherwise}
    \end{cases}
  \end{equation*}
\end{proposition}
\begin{proof}
  Let $x_0$ be the state where all nodes are disconnected
  and the internal state of their sites is set to $0$.
  We set $E(x_0) = 0$
  and construct a canonical path from $x_0$ to $x$
  to assign an energy $E(x)$ to a state $x$
  % by constructing a canonical path from $x_0$ to $x$
  as the sum of the energy contributions of each step along the path.
  The canonical path starts by performing all flips first
  (in the order set by $<_a$) and then all bindings.
  % In an FB-system there is at most
  % one transition between any two states.
  Using $E(y) - E(x) = \ln\,q(y,x) - \ln\,q(x,y)$ to compute
  the energy contribution of each flip and bind,
  we obtain \eqn{fb-energy}.
\end{proof}

Finally, we would like to know when an FB-system
has detailed balance.

\begin{proposition}
  An FB-system has detailed balance if and only if
  \eqns{flipflip}{flipbind} hold for all states,
  nodes and sites in the system.
\end{proposition}
\begin{proof}
  By \lems{flipflip}{flipbind},
  the forward direction has been already proved.
  The backward direction amounts to proving that
  the definition of \eqn{fb-energy}
  is independent of the choice of path,
  \ie that any alternative energy function $E'$
  that uses a non-canonical path for its definition
  is equivalent to $E$ as defined in \prop{fb-energy}.
  In a non-canonical path
  we might have two types of non-canonical flips:
  those that occur after a bind and those that happen in
  an order different to the one specified by $<_a$.
  The terms contributed by the latter can be transformed into
  canonical ones using \eqn{flipflip}
  and those contributed by the former using \eqn{flipbind}.
  % This last transformation is equivalent to swapping the order
  % in which the events occurred along the path.
  Non-canonical binds that occur before a flip
  can be transformed into canonical binds by a similar argument.
\end{proof}

% It is easy to see how an ANC model can be encoded as an FB-system
% by flattening the hierarchy of sites.

% Although a cooperative assembly system can be \emph{approximated}
% by an FB-system where the degree of a node is encoded in
% the number of sites in some fixed state, % in state 1,
% it is yet unclear if an FB-system can simulate
% a cooperative assembly system exactly.


%%% Local Variables:
%%% mode: latex
%%% TeX-master: "thesis"
%%% End: