Extreme robustness of scaling in sample space reducing processes explains Zipf's law in diffusion on directed networks

Think of the N possible states of a given system as stairs with different widths and imagine a ball bouncing downstairs with random step sizes. The probability of the downward bouncing ball to hit stair i is proportional to its width q(i), see figure 1(d). Given these prior probabilities q(i), the transition probability from stair j to stair i is

$$\begin{eqnarray}&&p(i| j)=\left\{\begin{array}{l}\tfrac{q(i)}{g(j-1)}\quad {\rm{if}}\ i\lt j\\ 0\ {\rm{otherwise}},\end{array}\right.\end{eqnarray} \tag{ 1 }$$

with $g(j-1)={\sum }_{{\ell }\lt j}q({\ell })$. Prior probabilities are normalised, ${\sum }_{i}q(i)=1$. We denote such a SSRP by ψ. One can safely assume the existence of a stationary visiting distribution, p, arising from many repetitions of process ψ and satisfying the following relation:

$$\begin{eqnarray}&&p(i)=\displaystyle \sum _{i\lt j\leqslant N}p(i| j)p(j).\end{eqnarray} \tag{ 2 }$$

Using equation (1), and forming the difference

$$\begin{eqnarray}&&\displaystyle \frac{p(i+1)}{q(i+1)}-\displaystyle \frac{p(i)}{q(i)}=-\displaystyle \frac{p(i+1)}{g(i)},\end{eqnarray} \tag{ 3 }$$

and by re-arranging terms we find that

$$\begin{eqnarray}&&\displaystyle \frac{p(i+1)g(i+1)}{q(i+1)}=\displaystyle \frac{p(i)g(i)}{q(i)},\end{eqnarray} \tag{ 4 }$$

where we use the fact that $g(i)+q(i+1)=g(i+1)$. Note that this is true for all values of i, and in particular

$$\begin{eqnarray}&&\displaystyle \frac{p(i)g(i)}{q(i)}=\displaystyle \frac{p(1)g(1)}{q(1)}=p(1),\end{eqnarray} \tag{ 5 }$$

since $g(1)=q(1)$. We arrive at the final result

$$\begin{eqnarray}&&p(i)=\displaystyle \frac{q(i)}{g(i)}p(1)\qquad {\rm{with}}\qquad \displaystyle \frac{1}{p(1)}=\displaystyle \sum _{j\leqslant N}\displaystyle \frac{q(j)}{g(j)}.\end{eqnarray} \tag{ 6 }$$

p(i) is the probability that we observe the ball ball bouncing downwards at stair i. Equation (6) shows that the path-dependence of the SSRP ψ deforms the prior probabilities of the states of a given system, $q(i)\to p(i)=\tfrac{q(i)}{g(i)}$. We can now discuss various concrete prior distributions. Note that equation (6) is exact and does not dependent on system size.

Polynomial priors and the ubiquity of Zipf's law: Given power law priors, $q(i)\sim {i}^{\alpha }$ with $\alpha \gt -1$, one can compute g up to a normalisation constant

$$\begin{eqnarray}&&g(i)=\displaystyle \sum _{j\leqslant i}{j}^{\alpha }=\displaystyle \frac{{i}^{\alpha +1}}{\alpha +1}+{ \mathcal O }({i}^{\alpha }),\end{eqnarray} \tag{ 7 }$$

which, when used in equation (6), asymptotically gives

$$\begin{eqnarray}&&p(i)\sim \displaystyle \frac{p(1)}{i},\end{eqnarray} \tag{ 8 }$$

i.e., Zipf's law. More generally, this result is true for polynomial priors, $q(j)\sim {\sum }_{i\leqslant m}{a}_{i}{j}^{\alpha (i)}$, where the degree of the polynomial $\alpha (m)=\max \{\alpha (i)\}$ is larger than −1, in the limit of large systems. Numerical simulations show perfect agreement with the theoretical prediction for various values of α, see figure 2(a) (circles, triangles, red squares).

Zoom In Zoom Out Reset image size

Fast decaying priors: The situation changes drastically for exponents $\alpha \lt -1$. For sufficiently fast decaying priors we have

$$\begin{eqnarray}&&g(i)\sim {\int }_{1}^{i}q(x){\rm{d}}x\sim g(1)=q(1).\end{eqnarray} \tag{ 9 }$$

The fast decay makes the contribution to g from large i's negligible. Under these circumstances equation (6) can be approximated for sufficiently large i's, as $p(i)\sim q(i)$. We encounter the remarkable situation that for fast decaying priors the SSRP, even though it is history dependent, follows the prior distribution. In this case the SSRP resembles a standard sampling process.

Exponential priors: For exponential priors, $q(i)\sim {{\rm{e}}}^{\beta {\rm{i}}}$, with $\beta \gt 0$, we find according to equation (6) that $p(i)=1/N$, i.e., a uniform distribution. To see this note that, up to a normalisation constant, g(i) is a geometric series

$$\begin{eqnarray*}&&g(i)=\displaystyle \sum _{j=1}^{i}{{\rm{e}}}^{\beta j}={{\rm{e}}}^{\beta }\displaystyle \frac{{{\rm{e}}}^{\beta {\rm{i}}}-1}{{{\rm{e}}}^{\beta }-1}.\end{eqnarray*}$$

Substituting it into equation (6), one finds the exact relation

$$\begin{eqnarray}&&p(i)=p(1)\displaystyle \frac{1-{{\rm{e}}}^{-\beta }}{1-{{\rm{e}}}^{-\beta {\rm{i}}}},\end{eqnarray} \tag{ 10 }$$

which can be safely approximated, for $i\gg 1$, by

$$\begin{eqnarray}&&p(i)\to p(1)\left(1-\displaystyle \frac{1}{{{\rm{e}}}^{\beta }}\right).\end{eqnarray} \tag{ 11 }$$

We observe that this is a constant independent of i. Accordingly, after normalisation, we will have $p(i)\sim 1/N$. Note that exponential priors describe a somewhat pathological situation. Given that a state i is occupied at time t, the probability to visit state $i-1$ is huge compared to all the other remaining states, so that practically all states will be sampled in a descending sequence: $i\to i-1\to i-2\to i-3\to \cdots 1$, which obviously leads to a uniform p. Again, numerical simulations show perfect agreement with the prediction, as shown in figure 2(a) (grey squares). Switching from polynomial to exponential priors, we switch the attractor from the Zipf's regime to the uniform distribution.

Noisy SSRPs are mixtures of a SSRP ψ and stochastic transitions between states that are not history-dependent. Following the previous scheme of the staircase picture, the noisy variant of the SSRP, denoted by ${\psi }_{\lambda }$, starts at N and jumps to any stair $i\lt N$, according to the prior probabilities q(i). At i the process now has two options: (i) with probability λ the process continues the SSRP and jumps to any $j\lt i$, or, (ii) with probability $1-\lambda$ jumps to any point $j\lt N$, following a standard process of sampling without memory. $1-\lambda$ is the noise strength. The process stops when stair 1 is hit. The transition probabilities for ${\psi }_{\lambda }$ read

$$\begin{eqnarray}&&p(i| j)=\left\{\begin{array}{l}\lambda \tfrac{q(i)}{g(j-1)}+(1-\lambda )q(i)\ \ {\rm{if}}\ i\lt j\\ (1-\lambda )q(i)\ {\rm{otherwise}}.\end{array}\right.\end{eqnarray} \tag{ 12 }$$

Note that the noise allows moves from j to i, even if $i\gt j$. Proceeding exactly as before we get

$$\begin{eqnarray}&&\displaystyle \frac{{p}_{\lambda }(i+1)}{q(i+1)}\left(1+\lambda \displaystyle \frac{q(i+1)}{g(i)}\right)=\displaystyle \frac{{p}_{\lambda }(i)}{q(i)},\end{eqnarray} \tag{ 13 }$$

where ${p}_{\lambda }(i)$ depicts the probability to visit state i in a noisy SSRP with parameter λ. As a consequence we obtain:

$$\begin{eqnarray}&&{p}_{\lambda }(i)={p}_{\lambda }(1)\displaystyle \frac{q(i)}{q(1)}\displaystyle \prod _{1\lt j\leqslant i}{\left(1+\lambda \displaystyle \frac{q(j)}{g(j-1)}\right)}^{-1}.\end{eqnarray} \tag{ 14 }$$

The product term can be safely approximated by

$$\begin{eqnarray}\displaystyle \prod _{1\lt j\leqslant i}{(\cdots )}^{-1} & = & \exp \left[-\displaystyle \sum _{1\lt j\leqslant i}\mathrm{log}\left(1+\lambda \displaystyle \frac{q(j)}{g(j-1)}\right)\right]\\ & \approx & \exp \left[-\displaystyle \sum _{1\lt j\leqslant i}\lambda \displaystyle \frac{q(j)}{g(j-1)}\right]\\ & \approx & \exp \left[-\lambda \mathrm{log}\left(\displaystyle \frac{g(i)}{q(1)}\right)\right]\\ & = & {\left(\displaystyle \frac{g(i)}{q(1)}\right)}^{-\lambda },\end{eqnarray} \tag{ 15 }$$

where we used $q(j)\sim {{\rm{d}}g/{\rm{d}}x| }_{j}$ and $\mathrm{log}(1+x)\sim x$ for small x, assuming that $x=\lambda \tfrac{q(j)}{g(j-1)}\ll 1$. Finally, we get

$$\begin{eqnarray}&&{p}_{\lambda }(i)\sim \displaystyle \frac{{p}_{\lambda }(1)}{q{(1)}^{1-\lambda }}\left(\displaystyle \frac{q(i)}{g{(i)}^{\lambda }}\right),\end{eqnarray} \tag{ 16 }$$

where ${p}_{\lambda }{(1)/q(1)}^{1-\lambda }$ acts as the normalisation constant. λ plays the role of a scaling exponent. For $\lambda \to 1$ (no noise), p_λ recovers the standard SSRP ψ of equation (1). For $\lambda =0$, we recover the case of standard random sampling, $p\to q$. It is worth noting that continuous SSRP display the same scaling behaviour (see appendix A). The particular case of $q(i)=1/N$ that was studied in [13], shows that λ turns out to be the scaling exponent of the distribution ${p}_{\lambda }(i)\sim 1/{i}^{\lambda }$. Note that these are not frequency- but rank distributions. They are related, however. The range of exponents $\lambda \in (0,1]$ in rank, represents the respective range of exponents $\alpha \in [2,\infty )$ in frequency, see e.g. [14] and appendix B. For polynomial priors, $q(i)\sim {i}^{\alpha }$ ($\alpha \gt -1$), one finds

$$\begin{eqnarray}&&{p}_{\lambda }(i)\sim {i}^{\alpha (1-\lambda )-\lambda }.\end{eqnarray} \tag{ 17 }$$

The excellent agreement of these predictions with numerical experiments is shown in figure 2(b). Finally, for exponential priors $q(i)\sim {{\rm{e}}}^{\beta {\rm{i}}}$ ($\beta \gt 0$) the visiting probability of for the noisy SSRP ${\psi }_{\lambda }$ becomes $p(i)\sim {{\rm{e}}}^{(1-\lambda )\beta {\rm{i}}}$, see table 1. Clearly, the presence of noise recovers the prior probabilities in a fuzzy way, depending on the noise levels. The following table sumarizes the various scenarios for the distribution functions p(i) for the different prior distributions q(i) and noise levels.

Table 1.  Distribution functions p(i) of SSRPs for the various prior distributions q(i). SSRP distributions with a noise level of $(1-\lambda )$ are indicated by ${p}_{\lambda }(i)$.

Prior	(sub-) logarithmic	Polynomial	Exponential
q(i)	${i}^{\alpha }$ ($\alpha \lt -1$)	${i}^{\alpha }$ ($\alpha \gt -1$)	${{\rm{e}}}^{\beta {\rm{i}}}$
p(i)	${i}^{\alpha }$	${i}^{-1}$	$\tfrac{1}{N}$
${p}_{\lambda }(i)$ noise	${i}^{\alpha }$	${i}^{\alpha (1-\lambda )-\lambda }$	${{\rm{e}}}^{(1-\lambda )\beta {\rm{i}}}$

A diffusion process on ${{ \mathcal G }}^{D}$ is now carried out by placing random walkers on the nodes randomly, and letting them take steps following the arrows in the network. They diffuse according to the weights in the network until they hit a target node and are then removed. We record the number of visits to all nodes and sort them according to the number of visits, obtaining a rank distribution of visits⁵ . We show the results from numerical experiments of 10⁷ random walkers on various DAGs in figure 4. In figures 4(a) and (b) we plot the rank distribution of visits to nodes for weighted Erdős–Rényi (ER) DAG networks. A weight w_ik is randomly assigned to each link ${e}_{{ik}}\in E$ from a given weight distribution p(w). Weights either follow a Poisson distribution, figure 4(a), or a power-law distribution, figure 4(b). In both cases Zipf's law is obtained in the rank distribution of node visits. For the same network we introduce noise with $\lambda =0.5$ and carry out the same diffusion experiment. The observed slope corresponds nicely with the predicted value of λ, as shown in figure 4(a) (red squares) for the Poisson weights.

Zoom In Zoom Out Reset image size

We computed rank distributions of node visits from diffusion on more general network topologies. In figure 4(c) we show the rank distribution of node visits where the substrate network is the citation network of high energy physics in the arXiv repository [33, 34], and the order is induced by the degree of nodes. Figure 4(d) shows the rank distribution of node visits from diffusion on an exponential DAG, that is generated by non-preferential attachment [35], where the order of nodes is again induced according to the degree. Both networks show Zipf's law in the rank distribution of node visits. This is remarkable since both networks are drastically different in topological terms.

For diffusion on random DAGs it is possible to obtain analytic results that are identical to equation (1), showing that Zipf's law is generally present in targeted diffusion.

We first focus on the definition of the prior probabilities in the context of diffusion on undirected networks. As stated above, q(i) is the probability that state i is visited in a random sampling process, see figures 1(b) and (c). In the network context this corresponds to the probability that node i is visited by a random walker. Assume that we have an undirected random graph ${ \mathcal G }(V,E)$ and that the N nodes are labelled $1,...N$. The probability that a random walker arrives at node i from a randomly chosen link of E, the network-prior probability of node i, is easily identified as

$$\begin{eqnarray}&&{q}_{G}(i)\equiv \displaystyle \frac{{k}_{i}}{2| E| },\end{eqnarray} \tag{ 18 }$$

where $| E|$ is the number of links in the graph; the factor 2 appears because a link contains 2 endpoints. If ${\sigma }_{G}\equiv \{{k}_{1},...,{k}_{N}\}$ denotes the undirected degree sequence q_G, is a simple rescaling of ${\sigma }_{G}$, i.e., ${q}_{G}=\tfrac{1}{2| E| }{\sigma }_{G}$. Using the same notation as before, the cumulative network-prior probability distribution is ${g}_{G}(i)\equiv {\sum }_{{\ell }\leqslant i}{q}_{G}({\ell })$.

From equation (18) and by assuming that in sparse graphs the probability of self-loops vanishes, i.e., $p({e}_{{ii}})\to 0$, one can compute the probability that a link e_ij exists in ${ \mathcal G }$, [32]

$$\begin{eqnarray}&&p({e}_{{ij}}\in E)=\displaystyle \frac{k(i)k(j)}{{\displaystyle \sum }_{{\ell }\leqslant N}k({\ell })}=2| E| {q}_{G}(i){q}_{G}(j),\end{eqnarray} \tag{ 19 }$$

where the second step is possible since ${\sum }_{{\ell }\leqslant N}k({\ell })=2| E|$. With this result, the out-degree of node labelled i in the graph ${{ \mathcal G }}^{D}$ can be approximated by

$$\begin{eqnarray}{k}_{i}^{{\rm{out}}} & = & \displaystyle \sum _{j\lt i}p({e}_{{ij}}\in E)\\ & = & 2| E| \displaystyle \sum _{j\lt i}{q}_{G}(i){q}_{G}(j)\\ & = & 2| E| {q}_{G}(i)\displaystyle \sum _{j\lt i}{q}_{G}(j)\\ & = & 2| E| {q}_{G}(i){g}_{G}(i-1).\end{eqnarray} \tag{ 20 }$$

Note that to compute ${k}_{i}^{{\rm{out}}}$ we only need take into account the (undirected) links which connect i to nodes with a lower label $j\lt i$, according to the labelling used for the DAG construction outlined above.

We can now compute the probability that a random walker jumps from node i to node j on the DAG ${{ \mathcal G }}^{D}$,

$$\begin{eqnarray}&&{p}_{G}(j| i)=\left\{\begin{array}{l}p(j| i,{e}_{{ij}}\in E)p({e}_{{ij}}\in E)\ \ {\rm{if}}\ i\gt j\\ 0\ \ {\rm{otherwise}}.\end{array}\right.\end{eqnarray} \tag{ 21 }$$

This is the network analogue of equation (1). Here $p(j| i,{e}_{{ij}}\in E)$ is the probability that the random walker jumps from i to j given that $i\gt j$ and the link e_ij exists in ${ \mathcal G }$. Clearly, this probability is

$$\begin{eqnarray}p(j| i,{e}_{{ij}}\in E) & = & \displaystyle \frac{1}{{k}_{i}^{{\rm{out}}}}\\ & = & {(2| E| {q}_{G}(i){g}_{G}(i-1))}^{-1},\end{eqnarray} \tag{ 22 }$$

Using equations (19) and (22) in (21) we get

$$\begin{eqnarray}&&{p}_{G}(j| i)=\left\{\begin{array}{l}\tfrac{{q}_{G}(j)}{{g}_{G}(i-1)};\ {\rm{if}}\ i\gt j\\ 0\ \ {\rm{otherwise}},\end{array}\right.\end{eqnarray} \tag{ 23 }$$

which has the same form as equation (1). Note that this expression only depends on q_G, i.e. the degrees of nodes in the undirected (!) graph ${ \mathcal G }$. The solution of equation (23) is obtained in exactly the same way as before for equation (1), and the node visiting probability of targeted diffusion on random DAGs is

$$\begin{eqnarray}&&{p}_{\text{}}(i)\propto \displaystyle \frac{{q}_{G}(i)}{{g}_{G}(i)},\end{eqnarray} \tag{ 24 }$$

which is the network analog of equation (6).

We finally show the results for a DAG that is based on an ER graph. For an ER graph, by definition, the probability for a link to exist is a constant $r\in (0,1]$, and $p({e}_{{ij}}\in E)=r$. Again we label all nodes by $1,...,N$ and build a DAG ${{ \mathcal G }}_{{ER}}^{D}$ as described above. It is not difficult to see that the out-degree of node i is ${k}^{{\rm{out}}}(i)=(i-1)r$, and, using this directly in equation (21), we get

$$\begin{eqnarray}&&{p}_{G}(j| i)=\left\{\begin{array}{l}\tfrac{1}{i-1}\ \ {\rm{if}}\ i\gt j\\ 0\ \ {\rm{otherwise}},\end{array}\right.\end{eqnarray} \tag{ 25 }$$

which is the standard equation for a SSRP with uniform prior probabilities q, [13]. This means that for the ER graph q_G(i) is a constant and ${g}_{G}(i)\sim i$. Using this in equation (24), we find that the node visiting probability is exactly Zipf's law, with respect to the ordering used to build the DAG

$$\begin{eqnarray}&&{p}_{\text{}}(i)\propto {i}^{-1}.\end{eqnarray} \tag{ 26 }$$

Note that this result is independent of r and, therefore, of the average degree of the graph.

With the example of the staircase in mind, we can describe a SSRP ψ over a continuous sampling space, see figure A1 . We start in the extreme of the interval, x = N, and we choose any point of Ω following the probability density q. Suppose we land in $x\lt N$. Then, at time t = 1 we choose at random some point $x^{\prime} \in {{\rm{\Omega }}}_{x}$ following a probability density proportional to q. We run the process until a point $z\in {{\rm{\Omega }}}_{1}$ is reached. Then the process stops. The SSRP ψ can be described by the transition probabilities between the elements of $x,y\in {\rm{\Omega }}$ such that $y\gt 1$ as follows

$$\begin{eqnarray}&&p(x| y)=\left\{\begin{array}{l}q(x)/g(y)\ \ {\rm{iff}}\ x\lt y\\ 0\ {\rm{otherwise}},\end{array}\right.\end{eqnarray} \tag{ A.2 }$$

where g(y) is the cumulative density distribution evaluated at point y,

$$\begin{eqnarray}&&g(y)={\int }_{{{\rm{\Omega }}}_{y}}q(x){\rm{d}}x={\int }_{1}^{y}q(x){\rm{d}}x+f(1).\end{eqnarray} \tag{ A.3 }$$

Zoom In Zoom Out Reset image size

We are interested in the probability density p which governs the frequency of visits along Ω after the sampling process ψ. To this end, we start with the following self-consistent relation for p,

$$\begin{eqnarray}&&p(x)={\int }_{x}^{N}p(x| y)p(y){\rm{d}}y.\end{eqnarray} \tag{ A.4 }$$

Recall that the integration limits ${\int }_{x}^{N}$ represent the fact that a particular state x can only be reached from a state $y\gt x$. By differentiating this integral equation we obtain:

$$\begin{eqnarray}&&\displaystyle \frac{{\rm{d}}p}{{\rm{d}}x}=\displaystyle \frac{{\rm{d}}}{{\rm{d}}x}\left({\int }_{x}^{N}p(x| y)p(y){\rm{d}}y\right).\end{eqnarray} \tag{ A.5 }$$

In agreement to equation (A.2), $p(x| y)=q(x)/g(y)$ if $y\gt 1$ and $y\gt x$. Equation (A.5) can be expanded using the Leibniz rule:

$$\begin{eqnarray}\displaystyle \frac{{\rm{d}}p(x)}{{\rm{d}}x} & = & {\displaystyle \int }_{x}^{N}\displaystyle \frac{{\rm{d}}p(x| y)}{{\rm{d}}x}p(y){\rm{d}}y-\displaystyle \frac{q(x)}{g(x)}p(x)\\ & = & \displaystyle \frac{1}{q(x)}\displaystyle \frac{{\rm{d}}q(x)}{{\rm{d}}x}{\displaystyle \int }_{x}^{N}\displaystyle \frac{q(x)}{g(y)}p(y){\rm{d}}y-\displaystyle \frac{q(x)}{g(x)}p(x)\\ & = & \displaystyle \frac{1}{q(x)}\displaystyle \frac{{\rm{d}}q(x)}{{\rm{d}}x}p(x)-\displaystyle \frac{q(x)}{g(x)}p(x).\end{eqnarray} \tag{ A.6 }$$

This leads to a differential equation governing the dynamics of SSRPs under arbitrary prior probabilities q,

$$\begin{eqnarray}&&\displaystyle \frac{{\rm{d}}p(x)}{{\rm{d}}x}=\left(\displaystyle \frac{1}{q(x)}\displaystyle \frac{{dq}(x)}{{\rm{d}}x}-\displaystyle \frac{q(x)}{g(x)}\right)p(x).\end{eqnarray} \tag{ A.7 }$$

The above equation can be easily integrated in the interval $(1,N]$. Observing that equation (A.7) can be rewritten as

$$\begin{eqnarray}&&\displaystyle \frac{{\rm{d}}p(x)}{p(x)}=\left[\displaystyle \frac{{\rm{d}}}{{\rm{d}}x}\mathrm{log}\left(\displaystyle \frac{q(x)}{g(x)}\right)\right]{\rm{d}}x.\end{eqnarray} \tag{ A.8 }$$

One finds:

$$\begin{eqnarray}&&\mathrm{log}p(x)=\mathrm{log}\left(\displaystyle \frac{q(x)}{g(x)}\right)+\kappa ,\end{eqnarray} \tag{ A.9 }$$

κ being an integration constant to be determined by normalisation. The above equation has as a general solution for points $x\in (1,N]$

$$\begin{eqnarray}&&p(x)=\displaystyle \frac{1}{Z}\displaystyle \frac{q(x)}{g(x)},\end{eqnarray} \tag{ A.10 }$$

where Z is the normalisation constant

$$\begin{eqnarray}&&Z={\int }_{0}^{N}\displaystyle \frac{q(y)}{g(y)}{\rm{d}}y.\end{eqnarray} \tag{ A.11 }$$

This demonstrates how the prior probabilities q are deformed when sampled through the SSRP ψ in the region $x\in (1,N]$. This is the analogous to equation (6) of the main text.

Suppose the interval ${\rm{\Omega }}=(0,N]$ and let us define a probability density q on Ω as in equation (A.1). The noisy SSRP ${\psi }_{\lambda }$ starts at x = N and jumps to any point in $x^{\prime} \in {\rm{\Omega }}$, according to the prior probabilities q. From $x^{\prime}$ the system has two options: (i) with probability λ the process jumps to any $x^{\prime\prime} \in {{\rm{\Omega }}}_{x^{\prime} }$, i.e., ${\psi }_{\lambda }$ continues the SSRP we described above or, (ii) with probability $1-\lambda$, ${\psi }_{\lambda }$ jumps to any point $x^{\prime\prime} \in {\rm{\Omega }}$, following a standard sampling process. The process stops when it jumps to a member of the sink set, namely to a $x\leqslant 1$. The transition probabilities now read $(\forall y\gt 1)$,

$$\begin{eqnarray}&&p(x| y)=\left\{\begin{array}{l}\lambda q(x)/g(y)+(1-\lambda )q(x)\ \ {\rm{iff}}\ x\lt y\\ (1-\lambda )q(x)\ {\rm{otherwise}},\end{array}\right.\end{eqnarray} \tag{ A.12 }$$

Note that the noise enables the process to move from y to x, in spite $x\gt y$. As we did in equation (A.4), we can find a consistency relation for the probability density p_λ of visiting a given point of Ω along a noisy SSRP

$$\begin{eqnarray}&&{p}_{\lambda }(x)=\lambda {\int }_{x}^{N}p(x| y)p(y){\rm{d}}y+(1-\lambda )q(x).\end{eqnarray} \tag{ A.13 }$$

If we take the derivative

$$\begin{eqnarray*}\displaystyle \frac{{\rm{d}}{p}_{\lambda }(x)}{{\rm{d}}x} & = & \lambda \displaystyle \frac{{\rm{d}}}{{\rm{d}}x}\left({\displaystyle \int }_{x}^{N}p(x| y){p}_{\lambda }(y){\rm{d}}y\right)+(1-\lambda )\displaystyle \frac{{\rm{d}}q(x)}{{\rm{d}}x}\\ & = & \lambda \displaystyle \frac{{\rm{d}}q(x)}{{\rm{d}}x}{\displaystyle \int }_{x}^{N}\displaystyle \frac{{p}_{\lambda }(y)}{g(y)}{\rm{d}}y-\lambda \displaystyle \frac{q(x)}{g(x)}{p}_{\lambda }(x)+(1-\lambda )\displaystyle \frac{{\rm{d}}q(x)}{{\rm{d}}x}\\ & = & \displaystyle \frac{\lambda }{q(x)}\displaystyle \frac{{\rm{d}}q}{{\rm{d}}x}{\displaystyle \int }_{x}^{N}\displaystyle \frac{q(x)}{g(y)}{p}_{\lambda }(y){\rm{d}}y-\lambda \displaystyle \frac{q(x)}{g(x)}{p}_{\lambda }(x)+(1-\lambda )\displaystyle \frac{{\rm{d}}q(x)}{{\rm{d}}x}\\ & = & \displaystyle \frac{1}{q(x)}\displaystyle \frac{{\rm{d}}q(x)}{{\rm{d}}x}({p}_{\lambda }(x)-(1-\lambda )q(x))-\lambda \displaystyle \frac{q(x)}{g(x)}{p}_{\lambda }(x)+(1-\lambda )\displaystyle \frac{{\rm{d}}q(x)}{{\rm{d}}x}\\ & = & \displaystyle \frac{1}{q(x)}\displaystyle \frac{{\rm{d}}q(x)}{{\rm{d}}x}{p}_{\lambda }(x)-\lambda \displaystyle \frac{q(x)}{g(x)}{p}_{\lambda }(x),\end{eqnarray*}$$

where the fourth step is performed taking the definition of ${p}_{\lambda }(x)$ given in equation (A.13). We therefore have the following differential equation for ${p}_{\lambda }(x)$,

$$\begin{eqnarray}&&\displaystyle \frac{{\rm{d}}{p}_{\lambda }(x)}{{\rm{d}}x}=\left(\displaystyle \frac{1}{q(x)}\displaystyle \frac{{\rm{d}}q(x)}{{\rm{d}}x}-\lambda \displaystyle \frac{q(x)}{g(x)}\right){p}_{\lambda }(x),\end{eqnarray} \tag{ A.14 }$$

which can be rewritten as

$$\begin{eqnarray*}&&\displaystyle \frac{{\rm{d}}{p}_{\lambda }(x)}{{p}_{\lambda }(x)}=\displaystyle \frac{{\rm{d}}}{{\rm{d}}x}\mathrm{log}\left(\displaystyle \frac{q(x)}{{g}^{\lambda }(x)}\right){\rm{d}}x.\end{eqnarray*}$$

Integrating it overall $x\in (1,N]$, we obtain

$$\begin{eqnarray}&&{p}_{\lambda }(x)=\displaystyle \frac{1}{{Z}_{\lambda }}\displaystyle \frac{q(x)}{{g}^{\lambda }(x)},\end{eqnarray} \tag{ A.15 }$$

which again demonstrates how the noisy SSRP deforms the underlying prior probabilities q, Z_λ being the normalisation constant. Interestingly, if $\lambda \lt 1$, i.e., if we consider a noisy SSRP, λ has the role of a scaling exponent. We observe that we recover the standard SSRP ψ described above in equation (A.2) if $\lambda \to 1$ (no noise) and the Bernouilli process following the prior probabilities q if we have total noise, as expected. The results for the continuous SSRPs are similar to the discrete case; compare equation (A.15) and equation (16).