Exact and Approximate Constants of Motion in Stochastic Contact Processes

1. Introduction

Constants of motion are central to the analysis of dynamical processes because they constrain the regions of state space that a system can access, thereby simplifying its mathematical description [1,2]. If a physical system admits one or more invariants –be they conserved energies, masses, charges, or more abstract integral quantities– its dynamics become confined to lower-dimensional manifolds within the original state space. This reduction is not merely computationally convenient, but provides structural insight into the qualitative behavior of the system. In deterministic dynamics, constants of motion can often be used to construct reduced models or to identify families of trajectories parameterized by invariant quantities. Even in systems that are not fully integrable, partial sets of invariants can help isolate slow manifolds, determine stability boundaries, and identify geometrical constraints governing long-term behavior [3]. Although mechanical systems supply the canonical examples of conservation laws –energy, linear and angular momentum– this specific framework is only one instance of a far broader methodological scaffold. Noether’s theorem [4], which links conserved mechanical quantities to symmetries, remains important as a conceptual foundation but, in many other contexts, constants of motion are not tied to continuous symmetries in any strict sense.

Beyond classical mechanics, constants of motion appear in a wide variety of domains where they play a similarly structural, though sometimes subtler, role. In stochastic processes, for example, conservation laws often take the form of invariant probability fluxes, conserved moments, or linear combinations of state variables that remain unchanged under the stochastic evolution [5,6]. These invariants can be crucial for deriving stationary distributions, identifying metastable sets, and reducing chemical reaction networks or population models to lower-dimensional effective descriptions. In chemical kinetics, conserved quantities restrict the dynamics of reaction networks to affine subspaces, allowing substantial dimensionality reductions even when the underlying reaction kinetics are nonlinear [7,8]. In the physics of fluids and kinetic theory, generalized Hamiltonian structures give rise to Casimir invariants, which constrain the admissible evolution of distribution functions and help classify equilibria and stability properties [9,10]. Across all these settings, constants of motion act as organizing principles: they delimit the possible changes of the system along time, restrict its long-term behavior, and offer a unifying language for reduction, stability analysis, and qualitative classification.

In this contribution, we focus on the role of constants of motion in the analysis of stochastic contact processes. These are aleatory processes defined on systems formed by agents, where the state of each agent changes with time as a result of the interaction with other peers [11]. As referenced below, contact processes play a fundamental role in the mathematical modeling of a wide class of phenomena in biological and socioeconomic systems, covering disease dissemination, information transmission, and commercial or financial transactions, among many others. Here, for a variety of stochastic contact processes related to the modeling of rumor and disease propagation, we construct linear constants of motion and show how they help obtaining information on the evolution of the system. We demonstrate that, while some processes admit exact nontrivial constants of motion, others can only have approximated conserved quantities, that fluctuate randomly around their initial value as time elapses. In all cases, we compare with the conserved quantities known to exist in the mean-field version of each process. We also explore the possibility of finding nonlinear constants of motion.

2. Exact Constants of Motion in Maki-Thompson Rumor Propagation

Stylized models of rumor propagation have become one of the main applications of contact processes in the mathematical approach to social phenomena. Starting from a handful of abstract stochastic processes proposed several decades ago [12,13,14], these models have been extended to encompass a wide class of factors that influence information transfer in social systems, from mass-media dissemination [15], to networked structures [16,17], to population heterogeneity [18].

Among the most traditional models, the Maki-Thompson process [13] emulates the propagation of a rumor over a homogeneous population formed by N agents. Evolution goes on by discrete steps. At any step, each agent can be in one of three states which, by analogy with epidemiological models (see Section 3), are named susceptible (S), infected (I), and refractory (or recovered, R). Susceptible agents have not heard the rumor yet, infected agents have heard the rumor and are willing to transmit it, and refractory agents have heard the rumor but have lost interest and do not transmit it. At each evolution step, one of the infected agents is selected at random and a second agent is chosen, also at random, from the remaining of the population. If the second agent is susceptible, it becomes infected, representing the transmission of the rumor. If, on the other hand, the second agent is infected or refractory, the first agent looses interest in the rumor and becomes refractory. Typically, the process is run starting from an initial condition with many susceptible agents (for instance, $S_0=N-1$), a few infected agents (for instance, $I_0=1$), and no refractory agents ($R_0=0$). Evolution ends when no more infected agents are left in the population. The first column of Table 1 shows the two possible events in the format of chemical reactions. In the second column, we quote the corresponding probabilities of each reaction at step n, where $S_n$ is the number of susceptible agents at that step.

2.1. Constants of Motion

For each event, we have constants of motion associated with each state S, I, and R. Since in a transmission event (first row of Table 1), a susceptible agent becomes infected, the number of susceptible agents, $S_n$, plus the index of the evolution step, n, remains constant. Likewise, the difference $I_n-n$ is a constant of motion. Finally, since the number of refractory agents does not change, $R_n$ remains constant as well. These three constants of motion are shown in the third column of the table. Analogous arguments show that, for an event where an infected agent becomes refractory (second row of Table 1), the associated constants of motion are $S_n$, $I_n+n$, and $R_n-n$. Naturally, moreover, any combination of the three constants of motion is also a constant of motion for each event.

Can we combine the constants of motion quoted above in such a way as to obtain a function of $S_n$, $I_n$, $R_n$, and n that is a conserved quantity for both events? To answer this question, we consider linear combinations of the constants of motion, requiring that the identity

$${\alpha }_1(S_n+n)+{\beta }_1(I_n-n)+{\gamma }_1{R}_n={\alpha }_2{S}_n+{\beta }_2(I_n+n)+{\gamma }_2(R_n-n)$$

(1)

holds for a suitable choice of the coefficients. In order to get a solution with constant coefficients, we must have ${\alpha }_1={\alpha }_2\equiv \alpha$, ${\beta }_1={\beta }_2\equiv \beta$, and ${\gamma }_1={\gamma }_2\equiv \gamma$. Equation (1) then implies

$$\alpha -2\beta +\gamma =0.$$

(2)

This linear equation has two independent solutions, for instance, $\alpha =\beta =\gamma =1$, and $\alpha =2$, $\beta =1$, $\gamma =0$. Replacing in the combinations of Eq. (1), we get

$$\begin{array}{ccc}\hfill K_1& =& S_n+I_n+R_n,\hfill \end{array}$$

(3)

$$\begin{array}{ccc}\hfill K_2& =& 2S_n+I_n+n,\hfill \end{array}$$

(4)

which are two quantities conserved in both events. In the first of these joint constants of motion, we recognize the total number of agents, $K_1=N$, which is trivially conserved due to the fact that no agent is created or destroyed in any of the events. On the other hand, $K_2$ is a nontrivial conserved quantity, combining the number of susceptible and infected agents, and the index of the evolution step [19]. Its value depends on the specific initial condition chosen to run the process: $K_2=2S_0+I_0$. In the following, we show how $K_2$ helps to get information on the evolution.

2.2. Solution for the Mean Number of Agents in Each State

Since the probability of each event (second column of Table 1) depends on $S_n$ only, it is possible to write a master equation for the probability that the number of susceptible agents is S as step n, $P_n\left(S\right)$, disregarding the number of infected and refractory agents. It reads

$$P_{n+1}\left(S\right)={\displaystyle \frac{S+1}{N-1}}{P}_n(S+1)+\left(1-{\displaystyle \frac{S}{N-1}}\right)P_n\left(S\right).$$

(5)

From this master equation, we immediately find the evolution law for the mean number of susceptible agents, ${\langle S\rangle }_n={\sum }_{S=0}^{N}S{P}_n\left(S\right)$,

$${\langle S\rangle }_{n+1}=\left(1-{\displaystyle \frac{1}{N-1}}\right){\langle S\rangle }_n,$$

(6)

whose solution reads

$${\langle S\rangle }_n={\left(1-{\displaystyle \frac{1}{N-1}}\right)}^n{S}_0.$$

(7)

Here, ${\langle ·\rangle }_n$ denotes the mean value, calculated at step n, over realizations of the stochastic process starting from the same initial condition $S_0$. Averaging Eq. () over realizations with the same initial condition and using Eq. (7), we get the mean number of infected agents:

$${\langle I\rangle }_n=I_0+2S_0\left[1-{\left(1-{\displaystyle \frac{1}{N-1}}\right)}^n\right]-n.$$

(8)

The mean number of recovered agents is just ${\langle R\rangle }_n=N-{\langle S\rangle }_n-{\langle I\rangle }_n$.

The solution for ${\langle I\rangle }_n$, Eq. (8), makes it possible to compute, although not always analytically, some important average features of the stochastic process. For instance, in the limit of large N, where ${[1-{(N-1)}^{-1}]}^n\approx exp(-n/N)$, we find the number of steps at which the number of infected agents reaches a maximum, and the value of ${\langle I\rangle }_n$ at the maximum, respectively,

$$n_{\max }=Nln{\displaystyle \frac{2S_0}{N}},{\langle I\rangle }_{\max }=I_0+2S_0-N\left(1+ln{\displaystyle \frac{2S_0}{N}}\right).$$

(9)

Note that these results make sense if $S_0>N/2$. Otherwise, the number of infected agents does not reach a maximum during the evolution. By equating ${\langle I\rangle }_n=0$, meanwhile, we can find implicit equations for the expected evolution step at which the process ends, and for the corresponding average number of susceptible agents. In the limit of large N, we respectively have

$$n_{\mathrm{end}}=I_0+2S_0\left[1-exp\left(-{\displaystyle \frac{n_{\mathrm{end}}}{N}}\right)\right],{\langle S\rangle }_{\mathrm{end}}=S_{0}exp\left(-{\displaystyle \frac{n_{\mathrm{end}}}{N}}\right).$$

(10)

For the initial condition $S_0\to N$ and $I_0\to 0$, we find $n_{\mathrm{end}}\approx 1.594N$, and obtain the classical result for “the proportion of the population never hearing a rumor” [14,19], ${\langle S\rangle }_{\mathrm{end}}=0.203N$.

2.3. Mean-Field Evolution

Mean-field differential equations for the Maki-Thompson model are found in the limits of a large number of agents and continuous time. Firstly, we note that, in continuous time, the sequential choice of infected agents of the stochastic model must be replaced by a scheme where, on the average, each infected agent has the opportunity of interacting with the rest of the population once per time unit. This is achieved if the duration of each evolution step is proportional to the inverse of the number of infected agents, concretely,

$$t_{n+1}-t_n={\displaystyle \frac{\tau }{I_n}},$$

(11)

where $t_n$ is the time of occurrence of step n, $I_n$ is the corresponding number of infected agents, and $\tau$ is a constant that fixes time units. Next, combining Eqs. (6) and (11), we estimate the time derivative of the number of susceptible agents as

$${\displaystyle \frac{\mathrm{d}S}{\mathrm{d}t}}\approx 〈{\displaystyle \frac{1}{t_{n+1}-t_n}}〉\left({\langle S\rangle }_{n+1}-{\langle S\rangle }_n\right)=-{\displaystyle \frac{1}{\tau N}}{\langle S\rangle }_n{\langle I\rangle }_n.$$

(12)

Within this approximation, we get

$$\begin{array}{c}\hfill {\displaystyle \frac{\mathrm{d}s}{\mathrm{d}t}}=-{\displaystyle \frac{si}{\tau }},\end{array}$$

(13)

where $s\left(t\right)={\langle S\rangle }_n/N$ and $i\left(t\right)={\langle I\rangle }_n/N$ are the average fractions of susceptible and infected agents, respectively. Proceeding in the same way for the time derivatives of the number of infected and refractory agents, we find

$${\displaystyle \frac{\mathrm{d}i}{\mathrm{d}t}}={\displaystyle \frac{si-i(i+r)}{\tau }},{\displaystyle \frac{\mathrm{d}r}{\mathrm{d}t}}={\displaystyle \frac{i(i+r)}{\tau }},$$

(14)

with $r\left(t\right)={\langle R\rangle }_n/N$.

Equations (13) and (14) provide the mean-field description of the Maki-Thompson model. Inspection immediately shows that the sum $s+i+r$ is a constant of motion. This represents, of course, the conservation of the total number of agents, thus corresponding to the constant of motion of Eq. (3). From the definition of the fractions s, i, and r, we have $s+i+r=1$. A second constant of motion can be obtained by calculating the ratio between the first of Eqs. (14) and Eq. (13):

$${\displaystyle \frac{\mathrm{d}i}{\mathrm{d}s}}=-{\displaystyle \frac{si-i(i+r)}{si}}=-2+{\displaystyle \frac{1}{s}},$$

(15)

where we have replaced $r=1-s-i$. Integrating with respect to s from reference values $s_0$ and $i_0$, we find the second constant for the mean-field equations:

$$K_{\mathrm{mf}}=2s+i-ln{\displaystyle \frac{s}{s_0}}\equiv 2s_0+i_0.$$

(16)

Remarkably, although the functional form of this conserved quantity does not coincide with that of $K_2$ in Eq. (), their numerical values in the mean-field limit are the same (up to a factor N): $K_2=N{K}_{\mathrm{mf}}$. In fact, identifying the reference values $s_0$ and $i_0$ with the initial condition, $s_0=S_0/N$ and $i_0=I_0/N$, we have

$$ln{\displaystyle \frac{s\left(t\right)}{s_0}}=ln{\displaystyle \frac{{\langle S\rangle }_n}{S_0}}=-{\displaystyle \frac{n}{N}},$$

(17)

where we have used Eq. (7) in the limit of large N. The two constants of motion of the mean-field equations play the same role as $K_1$ and $K_2$ in the Maki-Thompson model, and can be used to get the mean-field counterpart of the results found in Section 2.2 for the stochastic process.

2.4. Two Extensions of the Maki-Thompson Model

In this section, we analyze two extensions of the Maki-Thompson stochastic process, which complexify the model by adding more states accessible to agents, but can be treated along the same lines as in the preceding sections. The first extension introduces resistance to recovery for infected agents [19]. Specifically, an infected agent must attempt transmission to infected or refractory agents a total of M times before becoming in turn refractory. This variant is formulated by introducing M states of infection, ${\mathrm{I}}^{\left(m\right)}$ ($m=0,\dots ,M-1$), where m indicates the number of transmission events already attempted by the agents in that state. Table 2 summarizes the relevant events, probabilities, and constants of motion. The original model is recovered for $M=1$.

Proceeding as in Section 2.1, we can find M equations for the $M+2$ coefficients of a linear combination which gives a constant of motion for all the events of the model. We thus have two independent solutions. Not unexpectedly, one of them corresponds to the trivial constant of motion $N=S_n+{\sum }_{m=0}^{M-1}{I}_n^{\left(m\right)}+R_n$, namely, the total number of agents. The second, nontrivial joint constant of motion is

$$K=(M+1)S_n+\sum _{m=0}^{M-1}(M-m)I_n^{\left(m\right)}+n.$$

(18)

For this extension of the Maki-Thompson model, the probability of having S susceptible agents at step n, $P_n\left(S\right)$, also satisfies the master equation (5), so that the solution for ${\langle S\rangle }_n$ given by Eq. (7) still holds. Now, however, the existence of a constant of motion cannot be used to find the number of infected agents $I^{\left(m\right)}$ for each value of m, but rather the combination ${\sum }_{m=0}^{M-1}(M-m)I_n^{\left(m\right)}$.

The mean-field equations for the present model are

$$\begin{array}{cc}\hfill {\displaystyle \frac{\mathrm{d}s}{\mathrm{d}t}}& =-{\displaystyle \frac{s}{\tau }}\sum _{m=0}^{M-1}{i}^{\left(m\right)},\hfill \\ \hfill {{\displaystyle \frac{\mathrm{d}i}{\mathrm{d}t}}}^{\left(0\right)}& ={\displaystyle \frac{s}{\tau }}\sum _{m=0}^{M-1}{i}^{\left(m\right)}-{\displaystyle \frac{1-s}{\tau }}{i}^{\left(0\right)},\hfill \\ \hfill {{\displaystyle \frac{\mathrm{d}i}{\mathrm{d}t}}}^{\left(m\right)}& ={\displaystyle \frac{1-s}{\tau }}\left[i^{(m-1)}-i^{\left(m\right)}\right],m=1,\dots ,M-1,\hfill \\ \hfill {\displaystyle \frac{\mathrm{d}r}{\mathrm{d}t}}& ={\displaystyle \frac{1-s}{\tau }}{i}^{(M-1)},\hfill \end{array}$$

(19)

with $i^{\left(m\right)}=\langle I^{\left(m\right)}\rangle /N$, and s and r defined as in the preceding section. Direct calculation shows that, together with $s+{\sum }_{m=0}^{M-1}{i}^{\left(m\right)}+r=1$, the quantity

$$K_{\mathrm{mf}}=(M+1)s+\sum _{m=0}^{M-1}(M-m)i^{\left(m\right)}-ln{\displaystyle \frac{s}{s_0}}=(M+1)s_0+\sum _{m=0}^{M-1}(M-m)i_0^{\left(m\right)},$$

(20)

with $s_0$ and $i_0^{\left(m\right)}$ the initial values of the corresponding variables. For large N, in full analogy with the original Maki-Thompson model, the numerical value of $K_{\mathrm{mf}}$ coincides with that of K in Eq. (18) up to a factor N: $K=N{K}_{\mathrm{mf}}$.

In the second extension of the Maki-Thompson model considered here, which incorporates resistance to infection, a susceptible agent must be contacted M times by any infected agent to become infected. To formulate this variant, we introduce $S_n^{\left(m\right)}$ as the number of susceptible agents, that, at step n, have already been contacted m times by infected agents ($m=0,1,\dots ,M-1$). Events, probabilities, and constants of motion are listed in Table 3. The original model is reobtained for $M=1$.

As in the first extension of the Maki-Thompson process, two independent constants of motion can be found. The first, trivial one is the total number of agents, $N={\sum }_{m=0}^{M-1}{S}_n^{\left(m\right)}+I_n+R_n$. The second joint constant of motion reads

$$K=\sum _{m=0}^{M-1}(M-m+1)S_n^{\left(m\right)}+I_n+n.$$

(21)

Now, the master equations for the different sub-classes of susceptible agents are coupled to each other, and so are the equations for the corresponding average numbers of agents. Nevertheless, it is possible to find their general solution. For an initial condition where no susceptible agent has yet been contacted by an infected agent, $S_0^{\left(0\right)}=S_0$ and $S_0^{\left(m\right)}=0$ for $m=1,\dots ,M-1$, the solution reads

$${\langle S^{\left(m\right)}\rangle }_n=\left({\displaystyle \genfrac{}{}{0pt}{}{n}{m}}\right){\left(1-{\displaystyle \frac{1}{N-1}}\right)}^{n-m}{\displaystyle \frac{S_0}{{(N-1)}^m}}.$$

(22)

Replacing into Eq. (21) gives the average number of infected agents ${\langle I\rangle }_n$. With this result, it is possible to evaluate the step n at which the number of agents reaches a maximum, the height of the maximum, the total duration of the process, and the number of surviving susceptible agents, much as done in Section 2.2 for the original model. These calculations are not included here for the sake of brevity.

In the mean-field limit we have

$$\begin{array}{cc}\hfill {{\displaystyle \frac{\mathrm{d}s}{\mathrm{d}t}}}^{\left(0\right)}& =-{\displaystyle \frac{s^{\left(0\right)}i}{\tau }},\hfill \\ \hfill {{\displaystyle \frac{\mathrm{d}s}{\mathrm{d}t}}}^{\left(m\right)}& ={\displaystyle \frac{i}{\tau }}\left[s^{(m-1)}-s^{\left(m\right)}\right],m=1,\dots ,M-1,\hfill \\ \hfill {\displaystyle \frac{\mathrm{d}i}{\mathrm{d}t}}& ={\displaystyle \frac{s^{(M-1)}i}{\tau }}-{\displaystyle \frac{1}{\tau }}\left(1-\sum _{m=0}^{M-1}{s}^{\left(m\right)}\right)i,\hfill \\ \hfill {\displaystyle \frac{\mathrm{d}r}{\mathrm{d}t}}& ={\displaystyle \frac{1}{\tau }}\left(1-\sum _{m=0}^{M-1}{s}^{\left(m\right)}\right)i,\hfill \end{array}$$

(23)

with $s^{\left(m\right)}=\langle S^{\left(m\right)}\rangle /N$, and i and r defined as in the preceding section. The nontrivial constant of motion for these equations is

$$K_{\mathrm{mf}}=\sum _{m=0}^{M-1}(M-m+1)s_n^{\left(m\right)}+i_n-ln{\displaystyle \frac{s}{s_0}}=\sum _{m=0}^{M-1}(M-m+1)s_0^{\left(m\right)}+i_0.$$

(24)

In full analogy with the cases considered above, in the mean-field limit, the value of this constant of motion is related to that of Eq. (21) as $K=N{K}_{\mathrm{mf}}$.

3. Approximate Constants of Motion in SIR Epidemics

The SIR (susceptible-infected-recovered) model for epidemic propagation is among the most fundamental mathematical representations of the spreading of a contagious disease in a population of interacting agents. It belongs to a wide class of models where contact processes are used to describe the elementary events intervening in disease dissemination [20,21,22]. In the SIR model, the contact of an infected agent and a susceptible agent may lead the latter to become infected. In contrast with the Maki-Thompson model studied above, however, the transition from the infected state to the recovered state occurs spontaneously, without intervention of other agents. Recovered agents, in turn, cannot change their state.

Although the SIR model is frequently formulated in terms of mean field equations [22], it also admits several versions in the form of a stochastic contact process. By analogy with the Maki-Thompson process, at each time step, we here choose an infected agent at random. With probability u, the agent becomes recovered. With the complementary probability, $1-u$, another agent is chosen at random from the whole population and, if this second agent is susceptible, it becomes infected. Table 4 shows the elementary events of the model, their probabilities, and the corresponding constants of motion.

Proceeding as for the Maki-Thompson model to find joint conserved quantities for the stochastic process, we note that –since, now, the process is defined by three individual events– the coefficients in a linear combination of the individual constants of motion must satisfy two conditions,

$$\alpha -\beta =\beta -\gamma =0,$$

(25)

instead of the only condition of Eq. (2). Consequently, there is only one independent solution, $\alpha =\beta =\gamma =1$, corresponding to the trivial constant of motion $N=S_n+I_n+R_n$. Thus, the SIR stochastic process has no other conserved quantity.

We note, however, that the first event in Table 4, on the one hand, and the second and third events, on the other, have their own constants of motion, respectively,

$$K_1=S_n+I_n+n,K_2=S_n+I_n.$$

(26)

This suggests that a combination of $K_1$ and $K_2$, weighted with the total probabilities of the corresponding events,

$${\tilde{K}}_n=u{K}_1+(1-u)K_2=S_n+I_n+un,$$

(27)

may play a role similar to that of a constant of motion, although not being an exactly conserved quantity itself. To assess this conjecture, we explicitly calculate the stochastic evolution of ${\tilde{K}}_n$, finding

$${\tilde{K}}_{n+1}=\left\{\begin{array}{cc}{\tilde{K}}_n-(1-u),\hfill & \mathrm{with\; probability}u,\hfill \\ {\tilde{K}}_n+u,\hfill & \mathrm{with\; probability}1-u.\hfill \end{array}\right.$$

(28)

This stochastic process represents a one-dimensional symmetric random walk: its mean value and variance evolve as

$${\langle \tilde{K}\rangle }_{n+1}={\langle \tilde{K}\rangle }_n,{\sigma }_{n+1}^2={\sigma }_n^2+u(1-u),$$

(29)

respectively. In other words, ${\tilde{K}}_n$ is an approximate constant of motion, in the sense that its mean value is a conserved quantity. Its standard deviation, in turn, grows as $n^{1/2}$, like in any ordinary diffusion process. Note that, in combination with the trivial constant of motion N, ${\tilde{K}}_n$ straightforwardly yields the mean number of recovered agents over realizations with initial condition $R_0$:

$${\langle R\rangle }_n=R_0+un.$$

(30)

The approximate constant of motion ${\tilde{K}}_n$ acquires relevance when considering the mean-field version of the SIR model. To show this, first, we use the master equation of the process,

$$P_{n+1}(S,I)=u{P}_n(S,I+1)+(1-u){\displaystyle \frac{S+1}{N-1}}{P}_n(S+1,I-1)+\left[1-(1-u){\displaystyle \frac{S}{N-1}}\right]P_n(S,I),$$

(31)

to find evolution equations for the mean numbers ${\langle S\rangle }_n$ and ${\langle I\rangle }_n$:

$${\langle S\rangle }_{n+1}=\left(1-{\displaystyle \frac{1-u}{N-1}}\right){\langle S\rangle }_n,{\langle I\rangle }_{n+1}={\langle I\rangle }_n+(1-u){\displaystyle \frac{{\langle S\rangle }_n}{N-1}}-u.$$

(32)

The solution to the first of these equations is

$${\langle S\rangle }_n={\left(1-{\displaystyle \frac{1-u}{N-1}}\right)}^n{S}_0$$

(33)

for any fixed initial condition $S_0$. Next, from Eqs. (32) and operating as for the Maki-Thompson model, we get the mean-field description:

$${\displaystyle \frac{\mathrm{d}s}{\mathrm{d}t}}=-{\displaystyle \frac{(1-u)si}{\tau }},{\displaystyle \frac{\mathrm{d}i}{\mathrm{d}t}}={\displaystyle \frac{(1-u)si}{\tau }}-{\displaystyle \frac{ui}{\tau }},{\displaystyle \frac{\mathrm{d}r}{\mathrm{d}t}}={\displaystyle \frac{ui}{\tau }}.$$

(34)

Now, calculating $\mathrm{d}i/\mathrm{d}s$ and integrating with respect to s, we find a mean-field constant of motion:

$$K_{\mathrm{mf}}=s+i-{\displaystyle \frac{u}{1-u}}ln{\displaystyle \frac{s}{s_0}}=s_0+i_0.$$

(35)

Much as in Eq. (17), in the mean-field limit we can write

$$ln{\displaystyle \frac{s\left(t\right)}{s_0}}=ln{\displaystyle \frac{{\langle S\rangle }_n}{S_0}}=-{\displaystyle \frac{(1-u)n}{N}},$$

(36)

where we have used Eq. (33). Comparing Eq. (27), this shows that, in the limit, the approximate constant of motion ${\tilde{K}}_n$ coincides with the exact mean-field conserved quantity $K_{\mathrm{mf}}$ up to a multiplicative constant: ${\tilde{K}}_n=N{K}_{\mathrm{mf}}$.

The main panel of Figure 1 illustrates the evolution of the number of susceptible and infected agents along three realizations of the stochastic process, for a system with $N=100$ agents and recovery probability $u=0.3$, from an initial condition with one infected agent and no refractory agents. A comparison with the mean evolution given by Eqs. (32) is provided. The inset shows the evolution of the approximate constant of motion ${\tilde{K}}_n$, together with its expected mean value and standard deviation.

3.1. SIR Epidemics on Random Networks

An important extension of epidemiological models based on contact processes, aimed at approaching a more realistic representation of social structure, consists in considering populations distributed on a network [23,24,25]. In this variant, agents are situated at the network nodes, and network links represent their mutual contacts, so that infection can only occur between neighbor agents. The heterogeneity induced by the network of contacts on the population structure implies that susceptible and infected agents have more chances of, respectively, becoming infected and transmitting the infection when their number of neighbors –namely, their degree– is larger. The mathematical description of the process thus requires discerning between agents with different degrees. For each possible degree k, we denote by ${\mathrm{S}}^{\left(k\right)}$, ${\mathrm{I}}^{\left(k\right)}$, and ${\mathrm{R}}^{\left(k\right)}$, the corresponding number of agents in each state of the SIR epidemics.

The SIR stochastic contact process on networks is here implemented as follows. At each evolution step, we choose at random a network link connected to (at least) one infected agent. With probability u, this infected agent becomes recovered. With the complementary probability, $1-u$, if the second agent connected to the chosen link is susceptible, it becomes infected. Clearly, since one link is chosen at each step, the probability that any given infected or susceptible agent becomes involved in an event now depends on the respective degree k. The second column of Table 5 quotes the probabilities for each event in the case of a network with links distributed at random and with no statistical correlation between the degrees of neighbor agents.

Linear combinations of the constants of motion in the third column of Table 5 show that the only independent conserved quantities common to all the events are, trivially, the total number of agents with each degree k, $N^{\left(k\right)}=S_n^{\left(k\right)}+I_n^{\left(k\right)}+R_n^{\left(k\right)}$, for $k=1,2,\dots$. Meanwhile, a combination analogous to ${\tilde{K}}_n$ in Eq. (27), namely,

$${\tilde{K}}_n^{\left(0\right)}=S_n+I_n+un=\sum _k\left[S_n^{\left(k\right)}+I_n^{\left(k\right)}\right]+un,$$

(37)

where $S_n$ and $I_n$ are the total numbers of susceptible and infected agents at step n, satisfies the same stochastic evolution equation as ${\tilde{K}}_n$ in Eq. (28), and is thus an approximate constant of motion.

It turns out, however, that ${\tilde{K}}_n^{\left(0\right)}$ is just but one among infinitely many combinations that act as approximate constants of motion for the SIR model on networks. To see this, consider for example the quantities $Q_n^{\left(\gamma \right)}={\sum }_k{(k/z)}^{\gamma }[S_n^{\left(k\right)}+I_n^{\left(k\right)}]$, for any value of the exponent $\gamma$, where z is the mean number of neighbors per network site. In an event as in the first line of Table 5, where an infected agent becomes recovered at step n, we have $Q_{n+1}^{\left(\gamma \right)}=Q_n^{\left(\gamma \right)}-{(k_n/z)}^{\gamma }$, where $k_n$ is the degree of the infected agent involved in that event. On the other hand, for the events in the second and third rows, we have $Q_{n+1}^{\left(\gamma \right)}=Q_n^{\left(\gamma \right)}$. Averaging over the possible values of $k_n$ with their corresponding probabilities (second column of Table 5), it turns out that the quantity

$${\tilde{K}}_n^{\left(\gamma \right)}=Q_n^{\left(\gamma \right)}+u\sum _{m=0}^{n-1}{\mathrm{\Theta }}_m^{\left(\gamma \right)}=\sum _k{\left({\displaystyle \frac{k}{z}}\right)}^{\gamma }\left[S_n^{\left(k\right)}+I_n^{\left(k\right)}\right]+u\sum _{m=0}^{n-1}{\mathrm{\Theta }}_m^{\left(\gamma \right)},$$

(38)

with

$${\mathrm{\Theta }}_m^{\left(\gamma \right)}={\displaystyle \frac{{\sum }_k{k}^{\gamma +1}{I}_m^{\left(k\right)}}{z^{\gamma }{\sum }_{k}k{I}_m^{\left(k\right)}}},$$

(39)

is driven by the stochastic process

$${\tilde{K}}_{n+1}^{\left(\gamma \right)}=\left\{\begin{array}{cc}{\tilde{K}}_n^{\left(\gamma \right)}-(1-u){\mathrm{\Theta }}_n^{\left(\gamma \right)},\hfill & \mathrm{with\; probability}u,\hfill \\ {\tilde{K}}_n^{\left(\gamma \right)}+u{\mathrm{\Theta }}_n^{\left(\gamma \right)},\hfill & \mathrm{with\; probability}1-u;\hfill \end{array}\right.$$

(40)

cf. Eq. (28). For this process,

$${\langle {\tilde{K}}^{\left(\gamma \right)}\rangle }_{n+1}={\langle {\tilde{K}}^{\left(\gamma \right)}\rangle }_n,{\sigma }_{n+1}^2={\sigma }_n^2+u(1-u){\left[{\mathrm{\Theta }}_n^{\left(\gamma \right)}\right]}^2,$$

(41)

give the evolution of the mean value and the standard deviation of ${\tilde{K}}_n^{\left(\gamma \right)}$; cf. Eq. (29). In other words, ${\tilde{K}}_n^{\left(\gamma \right)}$ is an approximate constant of motion for any value of $\gamma$. Linear combinations of ${\tilde{K}}_n^{\left(\gamma \right)}$ for different values of $\gamma$ would produce infinitely many other approximate conserved quantities.

Mean-field equations for the stochastic SIR model on networks are obtained by assigning to each event a duration which appropriately takes into account the probability that an infected agent of degree k becomes involved in that event. Equation (11) now becomes

$$t_{n+1}-t_n={\displaystyle \frac{\tau k_n}{{\sum }_{k^{\prime }}{k}^{\prime }{I}_n^{\left(k^{\prime }\right)}}}\equiv {\displaystyle \frac{\tau k_n}{{\vartheta }_n}},$$

(42)

with ${\vartheta }_n$ defined as in the caption of Table 5. Here, $k_n$ is the degree of the infected agent involved at step n and, as in Eq. (11), $\tau$ fixes time units. The resulting equations for the fractions $s^{\left(k\right)}=\langle S^{\left(k\right)}\rangle /N$, $i^{\left(k\right)}=\langle I^{\left(k\right)}\rangle /N$, and $r^{\left(k\right)}=\langle R^{\left(k\right)}\rangle /N$ are

$$\begin{array}{cc}\hfill {{\displaystyle \frac{\mathrm{d}s}{\mathrm{d}t}}}^{\left(k\right)}& =-{\displaystyle \frac{(1-u)k{s}^{\left(k\right)}}{z\tau }}\sum _{k^{\prime }}{\displaystyle \frac{k^{\prime }}{z}}{i}^{\left(k^{\prime }\right)},\hfill \\ \hfill {{\displaystyle \frac{\mathrm{d}i}{\mathrm{d}t}}}^{\left(k\right)}& ={\displaystyle \frac{(1-u)k{s}^{\left(k\right)}}{z\tau }}\sum _{k^{\prime }}{\displaystyle \frac{k^{\prime }}{z}}{i}^{\left(k^{\prime }\right)}-{\displaystyle \frac{u{i}^{\left(k\right)}}{\tau }},\hfill \\ \hfill {{\displaystyle \frac{\mathrm{d}r}{\mathrm{d}t}}}^{\left(k\right)}& ={\displaystyle \frac{u{i}^{\left(k\right)}}{\tau }},\hfill \end{array}$$

(43)

for each degree k. Operating on these equations, it can be straightforwardly shown that the quantity

$$K_{\mathrm{mf}}=\sum _{k^{\prime }}{\displaystyle \frac{k^{\prime }}{z}}\left[s^{\left(k^{\prime }\right)}+i^{\left(k^{\prime }\right)}\right]-{\displaystyle \frac{u}{1-u}}{\displaystyle \frac{z}{k}}ln{\displaystyle \frac{s^{\left(k\right)}}{s_0^{\left(k\right)}}}$$

(44)

is a constant of motion, with $s_0^{\left(k\right)}$ the initial density of susceptible agents with k neighbors. We stress that, in spite of its functional form, the last term in the right-hand side of Eq. (44) is independent of k, as it becomes clear if the first of Eqs. (43) is rewritten as

$${\displaystyle \frac{1}{k}}{\displaystyle \frac{\mathrm{d}}{\mathrm{d}t}}ln{s}^{\left(k\right)}=-{\displaystyle \frac{1-u}{z\tau }}\sum _{k^{\prime }}{\displaystyle \frac{k^{\prime }{i}^{\left(k^{\prime }\right)}}{z}}.$$

(45)

Comparing Eqs. (44) and (38), we realize that the approximate constant of motion of the stochastic process to be related to $K_{\mathrm{mf}}$ is ${\tilde{K}}_n^{\left(1\right)}$. Up to a factor N, the summations in both constants are the same, so that the comparison must be carried out between the last term in the right-hand side of each equation. In the stochastic process, the mean evolution of susceptible agents is

$${\langle S^{\left(k\right)}\rangle }_n={\left[1-{\displaystyle \frac{(1-u)k}{z(N-1)}}\right]}^n{S}_0^{\left(k\right)}.$$

(46)

Replacing in the last term of Eq. (44), we get

$$-{\displaystyle \frac{u}{1-u}}{\displaystyle \frac{z}{k}}ln{\displaystyle \frac{s^{\left(k\right)}}{s_0^{\left(k\right)}}}={\displaystyle \frac{un}{N}},$$

(47)

in the large-N limit. As for the last term of Eq. (38) with $\gamma =1$, $u{\sum }_m{\mathrm{\Theta }}_m^{\left(1\right)}$, Eq. (39) implies that, if infected agents are evenly distributed over the network, we can approximate ${\mathrm{\Theta }}_m^{\left(1\right)}=\langle k^2\rangle /z^2$ for all m. Within this approximation, which should improve as the system size grows, we have

$$u\sum _{m=0}^{n-1}{\mathrm{\Theta }}_m^{\left(1\right)}=un{\displaystyle \frac{\langle k^2\rangle }{z^2}}.$$

(48)

Taking into account Eq. (46), it is clear that the expected proportionality ${\tilde{K}}_n^{\left(1\right)}=N{K}_{\mathrm{mf}}$ only holds if $\langle k^2\rangle =z^2$, i. e. for a regular random network [26], where all agents have the same number of neighbors. Note that this condition is verified, in particular, for a fully connected network, where $k=N-1$ for all agents. In this case, all agents can interact with each other, and we recover the SIR epidemic model without network.

Figure 2 shows, in the main panel, results for three realizations of the SIR stochastic process on Erdos-Rényi random networks [26] with $N=100$ agents, and recovery probability $u=0.2$, from an initial condition with one infected agent and no refractory agents. On the average, each agent has $z=10$ neighbors. Full curves stand for the evolution of the total number of susceptible and infected agents, respectively, $S_n$ and $I_n$. Dashed curves correspond to their expected values, averaged over realizations of both the evolution and the underlying network. The inset shows the evolution of the approximate constant of motion ${\tilde{K}}_n^{\left(1\right)}$, given by Eq. (38) for $\gamma =1$. The horizontal dashed line indicates its expected mean value, $zN$.

4. Nonlinear and “Failed” Constants of Motion

The conserved quantities derived for the stochastic processes discussed in the preceding sections were built as linear combinations of the constants of motion associated to each event, which we have quoted in Table 1 to Table 5. Nothing precludes, however, that a suitably chosen nonlinear combination of the individual constants of motion can fulfill that role for the whole process. To illustrate this situation, we consider a variant of the Maki-Thompson model of rumor propagation analyzed in Section 2, where the contact between two infected agents results in both agents becoming refractory. Any other event occurs as in the original model. Table 6 lists the events, their probabilities, and the individual constants of motion for this variant.

The only linear combination of the constants of motion in Table 6 which acts as a conserved quantity for the three events is, not unexpectedly, the total number of agents $N=S_n+I_n+R_n$. In principle, thus, we cannot rely on other combinations to get information on the evolution of the stochastic process. However, some algebra on the mean-field equations,

$${\displaystyle \frac{\mathrm{d}s}{\mathrm{d}t}}=-{\displaystyle \frac{si}{\tau }},{\displaystyle \frac{\mathrm{d}i}{\mathrm{d}t}}={\displaystyle \frac{si-i(2i+r)}{\tau }},{\displaystyle \frac{\mathrm{d}r}{\mathrm{d}t}}={\displaystyle \frac{i(2i+r)}{\tau }},$$

(49)

shows that

$$K_{\mathrm{mf}}={\displaystyle \frac{i+1}{s}}+2ln{\displaystyle \frac{s}{s_0}}\equiv {\displaystyle \frac{i_0+1}{s_0}}$$

(50)

is a constant of motion along the mean-field dynamics. This suggests that a constant of motion for the stochastic process, now involving a nonlinear combination of the number of agents in each state, may in fact exist.

To explore this possibility, we first write down the master equation for the probability $P_n(S,I)$ that, at step n, the numbers of susceptible and infected agents are S and I, respectively:

$$P_{n+1}(S,I)={\displaystyle \frac{S+1}{N-1}}{P}_n(S+1,I-1)+{\displaystyle \frac{I+1}{N-1}}{P}_n(S,I+2)+\left(1-{\displaystyle \frac{S+I}{N-1}}\right)P_n(S,I+1),$$

(51)

where each term in the right-hand side represents the contribution of one of the events in Table 6. From this master equation, we reobtain Eq. (6) for the evolution of ${\langle S\rangle }_n$. For ${\langle I\rangle }_n$, in turn, we get

$${\langle I\rangle }_{n+1}={\langle I\rangle }_n+{\displaystyle \frac{1}{N-1}}\left[{\langle S\rangle }_n-2({\langle I\rangle }_n-1)-{\langle R\rangle }_n\right],$$

(52)

with ${\langle R\rangle }_n=N-{\langle S\rangle }_n-{\langle I\rangle }_n$. Taking into account the evolution laws for ${\langle S\rangle }_n$ and ${\langle I\rangle }_n$, Eqs. (52) and (6), it turns out that the combination

$$\tilde{K}={\displaystyle \frac{{\langle I\rangle }_n+N-2}{{\langle S\rangle }_n}}-{\displaystyle \frac{2}{N-2}}n$$

(53)

is preserved along the evolution. Note however that the nature of this constant of motion is different from that of the examples considered in the preceding sections. In fact, it does not coincide with either a combination of the variables $S_n$ and $I_n$ or the mean value of any such combination. Rather, it is a nonlinear combination of mean values. In any case, in view of the solution for ${\langle S\rangle }_n$ given by Eq. (7), in the limit of large N, $\tilde{K}$ coincides with the mean-field constant of motion $K_{\mathrm{mf}}$ given by Eq. (50). Moreover, using the solution for ${\langle S\rangle }_n$ and the form of $\tilde{K}$, it is possible to find the solution for the mean number of infected agents, which reads

$${\langle I\rangle }_n=\left(I_0+N-2+{\displaystyle \frac{2n}{N-2}}{S}_0\right){\left(1-{\displaystyle \frac{1}{N-1}}\right)}^n-N+2.$$

(54)

Finally, we consider a complementary variant of the Maki-Thompson process, where an infected agent must necessarily contact a refractory agent to become in turn refractory. Namely, the contact between two infected agents produces no change. Table 7 shows the relevant parameters. As in the preceding variant, the only linear constant of motion is the total number of agents, $N=S_n+I_n+R_n$.

For this process, the mean-field equations read

$${\displaystyle \frac{\mathrm{d}s}{\mathrm{d}t}}=-{\displaystyle \frac{si}{\tau }},{\displaystyle \frac{\mathrm{d}i}{\mathrm{d}t}}={\displaystyle \frac{(s-r)i}{\tau }},{\displaystyle \frac{\mathrm{d}r}{\mathrm{d}t}}={\displaystyle \frac{ri}{\tau }}.$$

(55)

Taking the ratio of the first and third of these equations, it becomes apparent that the product

$$K_{\mathrm{mf}}=sr=s_0{r}_0$$

(56)

is a conserved quantity. In order to relate $K_{\mathrm{mf}}$ to a constant of motion for the stochastic process, it is convenient to use the master equation for the probability $P_n(S,R)$:

$$P_{n+1}(S,R)={\displaystyle \frac{S+1}{N-1}}{P}_n(S+1,R)+{\displaystyle \frac{R-1}{N-1}}{P}_n(S,R-1)+\left(1-{\displaystyle \frac{S+R}{N-1}}\right)P_n(S,R),$$

(57)

which yields

$${\langle SR\rangle }_{n+1}={\langle SR\rangle }_n.$$

(58)

In other words, the product ${\tilde{K}}_n=S_n{R}_n$ is an approximate conserved quantity for the stochastic process, whose mean value is proportional to the constant of motion of the mean-field equations: ${\langle \tilde{K}\rangle }_n=N^2{K}_{\mathrm{mf}}$.

As a final warning, we remark that the mean number of susceptible and refractory agents here evolve as

$${\langle S\rangle }_n={\left(1-{\displaystyle \frac{1}{N-1}}\right)}^n{S}_0,{\langle R\rangle }_n={\left(1+{\displaystyle \frac{1}{N-1}}\right)}^n{R}_0,$$

(59)

from initial conditions $S_0$ and $R_0$. If, operating as for the preceding variant in Eq. (53), we naively compute the product

$${\langle S\rangle }_n{\langle R\rangle }_n={\left[1-{\displaystyle \frac{1}{{(N-1)}^2}}\right]}^n{S}_0{R}_0,$$

(60)

we get a quantity that is not a constant of motion, but decreases monotonically as n grows. In other words, in the present case, a combination of stochastic variables inspired by a conserved quantity at the mean-field level, fails to give origin to a constant of motion. Only in the limit $N\to \infty$ does the procedure work, yielding a conserved quantity which, as expected, is proportional to that of the mean-field equations.

5. Conclusion

In this paper, we have considered a variety of stochastic contact processes where agents in a population change their individual state through random events of pair interaction. The focus has been put on the existence of conserved quantities, and on the role of these constants of motion in the mathematical analysis of such processes. Our procedure can be summarized as follows. For each elementary interaction event in the process, we first identify its individual constants of motion, which immediately derive from the changes in the state of the involved agents. Then, we find combinations of the individual constants of motion that are conserved quantities of all the interaction events and, therefore, are strictly conserved along the whole dynamics. In principle, these combinations could be arbitrary functions of the individual constants of motion, but we have mainly considered linear combinations, for which the procedure can be carried out systematically. Individual constants of motion can also be weighted by the probabilities of the corresponding events, and then linearly combined to get approximately conserved quantities, which fluctuate randomly around a well defined value as the system evolves. Exact and approximate constants of motion can finally be used to ease the analytical treatment by reducing the number of relevant independent variables. Comparison with their counterparts in a mean-field description of the dynamics also proves to be meaningful: in the appropriate limit, the values of stochastic and mean-field conserved quantities usually coincide with each other.

For the Maki-Thompson model of rumor propagation, which is defined by two possible interaction events with three accessible states, the above procedure yields two exact constants of motion. The first, trivial one is just the total number of agents. Together with the second constant of motion, it reduces the state space of the system to one single dimension which, using the master equation of the stochastic process, yields the full solution to the problem. Moreover, in the continuous-time limit for large systems, the values of the two constants of motion coincide with those found for the corresponding mean-field equations, even when their functional forms are not the same. Extensions of the Maki-Thompson model, adding intermediate states to the evolution of each agent, also contain two exact constants of motion. However, since the state space in now higher-dimensional, the usefulness of the conserved quantities in solving the problem is more limited.

In the SIR model of disease spreading, where agents can access three different states as well, the number of possible elementary events grows to three and, consequently, there is only one exact constant of motion, namely, the number of agents. However, it is still possible to find an approximate nontrivial conserved quantity. Its dynamics coincide with those of a symmetric random walker, so that its mean value remains constant while its standard deviation grows as the square root of the number of steps. Thus, this approximate constant of motion plays the same role, with respect to the evolution averaged over realizations, as the exact nontrivial conserved quantity of the Maki-Thompson model, and allows for the full solution to the stochastic problem in terms of the mean values of the relevant variables. Moreover, as in the case of rumor propagation, its value coincides with a constant of motion of the mean-field equations in the appropriate limit. Similar results are obtained for the SIR model evolving on random networks, although in this case the exact coincidence with the mean-field constant of motion is restricted to the case of regular networks.

Finally, we have considered a couple of stochastic processes in which the corresponding mean-field equations suggest the existence of conserved quantities given by nonlinear combinations of the constants of motion of the elementary events. In these cases, it was possible to show that a suitable extension of the mean-field conserved quantities to the level of the stochastic description does yield approximate constants of motion, although not all extensions successfully connect conserved quantities at the two levels.

Although the models for rumor and disease propagation considered in this paper provide a fair illustration of how conserved quantities in stochastic process yield useful information on their evolution, and give insight on the connection with mean-field descriptions, it is clear that these examples do not exhaust the possible uses of such techniques. Stochastic models for chemical reactions [5,27], biological population dynamics [28,29,30,31], and socioeconomic phenomena [11,32] are foreseeable fields of application of the same approach.