# Everlasting impact of initial perturbations on first-passage times of non-
Markovian random walks

N. Levernier1, T. V. Mendes2, O. Bénichou3, R. Voituriez3,4, T. Guérin2 1 Aix
Marseille Univ., Université de Toulon, CNRS, CPT, Turing Center for Living
Systems, 13009 Marseille, France  2Laboratoire Ondes et Matière d’Aquitaine,
University of Bordeaux, Unité Mixte de Recherche 5798, CNRS, F-33400 Talence,
France  3Laboratoire de Physique Théorique de la Matière Condensée, CNRS/UPMC,
4 Place Jussieu, 75005 Paris, France  4Laboratoire Jean Perrin, CNRS/UPMC, 4
Place Jussieu, 75005 Paris, France

###### Abstract

Persistence, defined as the probability that a fluctuating signal has not
reached a threshold up to a given observation time, plays a crucial role in
the theory of random processesMajumdar1999 ; ReviewBray . It quantifies the
kinetics of processes as varied as phase ordering, reaction diffusion or
interface relaxation dynamics. The fact that persistence can decay
algebraically with time with non trivial exponents has triggered a number of
experimental wong2001measurement ; tam1997first ; dougherty2002experimental ;
merikoski2003temporal ; soriano2009universality  and theoretical
Majumdar1996uq ; Derrida1996 ; derrida1995exact ; Majumdar1996 ;
Majumdar2001fk ; majumdar1996global ; ehrhardt2002series ; lee1997persistence
; constantin2003infinite ; Krug1997 ; Molchan1999  studies. However, general
analytical methods to calculate persistence exponents cannot be applied to the
ubiquitous case of non-Markovian systems relaxing transiently after an imposed
initial perturbation. Here, we introduce a theoretical framework that enables
the non perturbative determination of persistence exponents of $d$-dimensional
Gaussian non-Markovian processes with general non stationary dynamics relaxing
to a steady state after an initial perturbation. Two prototypical classes of
situations are analyzed: either the system is subjected to a temperature
quench at initial time, or its past trajectory is assumed to have been
observed and thus known. Altogether, our results reveal and quantify, on the
basis of Gaussian processes, the deep impact of initial perturbations on
first-passage statistics of non-Markovian processes. Our theory covers the
case of spatial dimension higher than one, opening the way to characterize
non-trivial reaction kinetics for complex systems with non-equilibrium initial
conditions.

The persistence $S(t)$ is the probability that a random process $x(t)$ has not
reached a threshold up to time $t$Majumdar1999 ; ReviewBray . This quantity is
a natural tool in non equilibrium statistical physics to probe the history of
various systems undergoing phase orderingderrida1994non ; derrida1995exact ;
Majumdar1996  or reaction diffusion dynamicsReviewBray , or to quantify the
efficiency of target search problems Redner:2001a ; Condamin2007 ;
Benichou2010 ; Schuss2007 ; godec2016universal ; delorme2015maximum ;
guerin2016mean ; levernier2019survival ; dolgushev2015contact ; Guerin2012a .
It has been recognized that the long time decay of persistence is often
algebraic, $S(t)\sim t^{-\theta}$, where the persistence exponent $\theta$ is
non trivial as soon as the process is non-Markovian (i.e. displays memory
effects).

As a matter of fact, even for seemingly simple Gaussian dynamics where all
correlation functions are known, $\theta$ is generally non-trivial and not
known in closed form. This has triggered an intense theoretical activity for
its determination. Existing approaches to quantify persistence exponents of
Gaussian processes can be classified according to the nature, stationary or
not, of the increments $x(t+\tau)-x(t)$. If these increments are stationary at
all times, meaning that their statistics do not depend on the observation time
$t$ (such as in the case of the fractional Brownian motion), $\theta$ is
exactly known Krug1997 ; Molchan1999 ; levernier2018universal  111In $d$
dimensions, for scale-invariant processes with stationary increments
levernier2018universal , $\theta=1-Hd$. In the opposite case where the
increments always depend on the observation time $t$ and thus never reach a
stationary dynamics (i.e. are stationary at no times), persistence exponents
have been calculated for the specific cases of the random acceleration process
BURKHARDT1993 ; deSmedt2001partial ; poplavskyi2018exact  or systems in which
the dynamics occurs at zero temperature poplavskyi2018exact ;
dornic2018universal ; derrida1994non ; derrida1995exact ; Majumdar1996 ;
Majumdar1996uq ; Derrida1996 ; watson1996persistence ; newman2001critical ,
and is thus deterministic with random initial conditions.

However, numerous physical situations display a relaxation dynamics –
typically after an initial perturbation – that becomes stationary only after a
transient regime. This is the rule for processes interacting with many degrees
of freedom, subjected to thermal fluctuations during the dynamics, but
prepared in a non-equilibrium or perturbed state. As a prototypical example,
consider a tagged monomer of a flexible polymer initially equilibrated at a
temperature $T\not=1$, and quenched to a different temperature $T_{0}=1$ at
time $t\geq 0$. The dynamics of the tagged monomer keeps transiently track of
this initial perturbation, and relaxes to the equilibrium state at $T_{0}$
with stationary increments. Persistence properties for such process with non
stationary increments (ie displaying aging), which in this example are
instrumental to quantify the reaction kinetics of the polymer with a given
reactive site, remain largely unknown. In fact, there is a fundamental reason
why standard methods to calculate persistence exponents cannot be applied for
transiently aging processes (see SI, sections A and B)222More precisely, for
transiently aging processes, the independent interval approximation, which is
usually applied to calculate the statistics of zero crossing of the Gaussian
process obtained after Lamperti transform, cannot be applied since intervals
between zero crossing become ill-defined, see SI Section A and B. . The only
available results for similar problems are limited to one-dimensional
processes and provide bounds for the persistence exponents as well as
perturbative expansions for weakly non-Markovian processesKrug1997 .

Here, we develop a general theoretical framework that enables the
determination of the persistence exponents of general Gaussian processes
displaying such transient aging dynamics. We stress that these Gaussian
processes are non-Markovian (display memory effects), and appear in a wide
range of contexts Min2005 ; kou2004generalized ; wei2000single ; turiv2013 ;
ochab2011scale ; cutland1995stock ; burnecki2012universal ;
ernst2012fractional ; weiss2013single ; mason1995optical . Our method enables
us to reveal and quantify the impact of initial conditions, such as a
temperature quench, on the persistence exponent. We also consider the case
where the past trajectory of the random walker is known, e.g. because it has
been observed. We show that the very observation of this past trajectory
modifies the persistence exponent which is quantified by our approach.
Importantly, our theory covers the physically relevant and widely unexplored
case of persistence for non-Markovian random walkers living in a space of
dimension higher than one333Note that here the spatial dimension $d$ is
different from the spatial dimension $D$ appearing in the problem of diffusive
persistenceMajumdar1996uq ; Derrida1996 , where the process studied is
$x(t)=\phi(\mathbf{0},t)$ where $\partial_{t}\phi=\nabla^{2}\phi$, $\nabla$
being the $D$-dimensional nabla operator. Here, higher dimensions means that
the variable $x$ it-self evolves in a $d$ dimensional space..

We first consider a one-dimensional isotropic non-Markovian Gaussian
stochastic process $x(t)$, which represents the position of a random walker at
time $t$. It is entirely defined by its mean value, assumed for simplicity to
be constant with time (unbiased process), and its covariance ${\rm
Cov}(x(t),x(t^{\prime}))=\sigma_{0}(t,t^{\prime})$. This covariance is assumed
to be given and to take the standard self-similar scalingReviewBray  form at
long times $t,t^{\prime}\gg 1$, $\sigma_{0}(t,t^{\prime})\sim
t^{2H}G(t/t^{\prime})\equiv\sigma(t,t^{\prime})$, where $H$ is the usual Hurst
exponent. We chose our units of time so that $G(1)=1$. At long times, the mean
square displacement $\sigma(t,t)=t^{2H}$ is assumed to diverge so that the
particle does not remain close to its initial position, which leads to $H>0$.
Furthermore, we assume that the statistics of the increments $x(t+\tau)-x(t)$
become stationary at long times, i.e. become independent of the observation
time $t$ when $t\to\infty$. This implies the existence of a transient regime
associated to the progressive decay of the memory of the initial state, and
defines a stationary covariance $\sigma_{s}$ given by

$\displaystyle\sigma_{s}(\tau,\tau^{\prime})=$
$\displaystyle\lim_{t\to\infty}\langle[x(t+\tau)-x(t)][x(t+\tau^{\prime})-x(t)]\rangle.$
(1)

Of note, the persistence exponent $\theta$ is known to be given by
$\theta=1-H$ under the stronger hypothesis that the statistics of the
increments is stationary at any time (i.e. when $\sigma_{s}=\sigma_{0}$)
Krug1997 ; Molchan1999 . The class of random walks that we consider here
covers a broad spectrum of non-Markovian processes used in physics, and in
particular both subdiffusive ($H<1/2$) and superdiffusive ($H>1/2$) walks.

Theoretical method to determine $\theta$.  Our starting point is the following
generalization of the renewal equation Redner:2001a

$\displaystyle p(0,t)=\int_{0}^{t}d\tau F(\tau)p(0,t|\mathrm{FPT}=\tau),$ (2)

which results from a partition over the first-passage event. In this equation,
$p(0,t)$ stands for the probability density that the random walker is at
position $x=0$ at time $t$, $F$ is the first-passage time (FPT) density and
$p(0,t|\mathrm{FPT}=\tau)$ is the probability density that $x=0$ at time $t$
given that the first-passage event occurred at time $\tau$.

To proceed further, we assume that the stochastic process in the future of the
FPT, defined by $y(t)\equiv x(t+\mathrm{FPT})$, is Gaussian with so far
undetermined mean $\mu(t)$ and covariance $\sigma_{\pi}(t,t^{\prime})$. Such
Gaussian approximation has proved successful to seize memory effects to
predict mean first-passage times of Gaussian random walkers in confinement
with stationary increments guerin2016mean ; levernier2019survival ;
levernier2020kinetics ; in the present context simulations show the broad
validity of this hypothesis (see SI, Fig. S1). A first result of our approach
is that the exponent $\theta$ is linked to the large time behaviour of
$\sigma_{\pi}(t,t)$, which is found from Eq. (2) to behave like (see SI,
Section D)

$\sigma_{\pi}(t,t)-\sigma_{0}(t,t)\underset{t\to\infty}{\propto}\
t^{2H-\theta}$ (3)

This means that the calculation of the exponent $\theta$ amounts to that of
the covariance $\sigma_{\pi}(t,t^{\prime})$ of the trajectories in the late
future of the first-passage.

Relying on a generalization of Eq. (2) to link the two-time joint probability
distribution functions of $x(t_{1}),x(t_{2})$ and the FPT density, we obtain a
self-consistent equation for the distribution of trajectories in the future of
the FPT, leading in the large time limit to (see SI, Section D for details):

$\displaystyle\int_{0}^{\infty}\frac{dt}{t^{H}}\Big{\\{}$
$\displaystyle\rho(t+\tau,t+\tau^{\prime})-\rho(t+\tau,t)\frac{\sigma(t+\tau^{\prime},t)}{\sigma(t,t)}-\rho(t+\tau^{\prime},t)\frac{\sigma(t+\tau,t)}{\sigma(t,t)}$
$\displaystyle+3\rho(t,t)\frac{\sigma(t+\tau,t)\sigma(t+\tau^{\prime},t)}{2\sigma(t,t)^{2}}-\frac{\rho(t,t)}{2\sigma(t,t)}\left[\sigma(t+\tau,t+\tau^{\prime})-\sigma_{K}(\tau,\tau^{\prime})\right]\Big{\\}}=0$
(4)

where $\rho\equiv\sigma_{\pi}-\sigma_{0}$ (for large times). Here,

$\displaystyle\sigma_{K}(t,t^{\prime})=\begin{cases}\sigma(t,t^{\prime})&\text{if
}\theta>1-H\\\ \sigma_{s}(t,t^{\prime})&\text{if }\theta<1-H\end{cases}$ (5)

Next, we find that the linear equation (Everlasting impact of initial
perturbations on first-passage times of non-Markovian random walks) admits
solutions of the scaling form
$\rho(t,t^{\prime})=t^{2H-\theta}z_{\theta}(t/t^{\prime})$, where
$z_{\theta}(u)$ satisfies a linear integral equation of the form

$\int_{0}^{1}K_{\theta}(u,v)\left[z_{\theta}(u)-z_{\theta}(1)\left(1-\frac{(2H-\theta)(1-u)}{2}\right)\right]{\rm
d}u=f_{\theta}(v),$ (6)

where $K_{\theta}$ and $f_{\theta}$ are given in SI (Section D) in terms of
$\sigma$. It is found that generic solutions $z_{\theta}(u)$ display
divergences for small $u$, and we argue that $\theta$ is obtained by imposing
that $z_{\theta}(u)$ is regular. We expect that this selection criterium is
valid at least for $2H-\theta>0$ since it amounts in this case to impose that
$\rho(t,0)=0$. Self-consistency reasons also lead us to restrict the analysis
to $H>1/3$ (see SI, Section E). In practice, the linear integral equation Eq.
(6) is solved numerically for a test value $\theta_{\rm test}$ and yields a
diverging solution $z_{\theta_{\rm test}}(u)\sim A(\theta_{\rm
test})u^{-\alpha(\theta_{\rm test})}$; the persistence exponent $\theta$ is
then obtained iteratively by enforcing that the prefactor vanishes, $A=0$ (see
SI, Section D). This finally provides a constructive, non perturbative
determination of the persistence exponent $\theta$ for Gaussian process with
general non stationary dynamics, which is the central result of this paper.

Applications.  We now show how these results enable us to determine the impact
of initial conditions in two physically relevant cases. The first type of
problems (called type I here) is the determination of $\theta$ in systems
which relax after a sharp temperature quench that occurs at initial time,
which is a very generic situation. Typically, physical realizations of the
random process $x(t)$ can be the position of a monomer in various models of
macromolecules or the local height of an interface, which span a number of
values of $H$. In all these models, assuming that the initial state for $t\leq
0$ is an equilibrium state at temperature $T\neq 1$, while the dynamics at
$t>0$ occurs at temperature $T_{0}=1$, the covariance function for
$t,t^{\prime}>0$ takes the form (see Ref.Krug1997  and SI, Section C)

$\displaystyle\sigma(t,t^{\prime})\propto T(t^{2H}+t^{\prime
2H})+(1-T)(t+t^{\prime})^{2H}-|t-t^{\prime}|^{2H}.$ (7)

Of note, the temperature $T$ before the quench can be lower or larger than the
temperature $T_{0}=1$ of the dynamics for $t>0$. Examples of survival
probabilities obtained from simulations are displayed in Fig 1(a), which
clearly shows that the persistence exponent depends on the choice of initial
conditions, and that the dependence of the persistence exponent on temperature
is correctly predicted by our approach [Fig 1(b)]. Remarkably, the values of
$\theta$ for different temperatures span a large set of values and are
markedly different from their value $\theta=1-H$ in the stationary state. In
one example of simulations, we also recorded the trajectories in the future of
the first-passage and measured numerically the function $z_{\theta}(x)$, which
shows good agreement with the theoretical prediction [see Fig 1(c)]. This
figure also illustrates our procedure to determine $\theta$ as defined above:
the calculated $z_{\theta}(u)$ show divergence for small $u$ whenever $\theta$
is above or below its exact value. In Fig. 1(d), we check that our theory is
also correct for different $H$ (focusing on $T=0$), for both superdiffusive
and subdiffusive processes, and even far from the Markovian regime $H=1/2$. In
addition, explicit results can be obtained by analyzing our formalism [Eq.
(Everlasting impact of initial perturbations on first-passage times of non-
Markovian random walks)] perturbatively in the limit $\varepsilon=H-1/2\to 0$.
An expansion up to second order leads for any temperature $T$ before the
quench to

$\displaystyle\theta_{\mathrm{I}}=1-H-2(\sqrt{2}-1)(1-T)(H-1/2)+\\{a_{1}\
(1-T)[a_{2}+(1-T)]\\}(H-1/2)^{2}+\mathcal{O}\left((H-1/2)^{3}\right),$ (8)

where analytical expressions of $a_{1}$ and $a_{2}$ are given in SI (Section
F), with numerical estimates $a_{1}\simeq 1.77$, $a_{2}\simeq 1.28$.
Interestingly, in the particular cases $T=0$ and $T=1$, the first order terms
coincide with the exact first order solution of Ref. Krug1997 , which points
towards the exactness of our approach at this order. These perturbative
results are in good agreement with simulation results [Fig. 1(c),(d)].
Finally, these results show that an imposed initial perturbation – here a
temperature quench, deeply impacts the first-passage statistics of the system.
In the case of subdiffusive (or antipersistent) dynamics ($H<1/2$, realized
typically in polymer models), it is found that, because of long range memory
effects, an initial quench from a high ($T>1$) to a low ($T_{0}=1$)
temperature can strongly slow down the first-passage kinetics ($\theta<1-H$),
while a quench from low to high temperatures accelerates the kinetics
($\theta>1-H$); opposite conclusions are reached for superdiffusive (or
persistent) dynamics ($H>1/2$).

Figure 1: Persistence for quenched fBM (type I problem) in $d=1$. (a) Example
of survival probabilities for interface dynamics with $H=3/8$ at different
temperatures. The slopes of the continuous lines is the value of $\theta$
predicted in our approach. (b) Systematic comparison of $\theta$ as measured
in simulations versus theoretical values for different $T$ and $H=3/8$. (c)
Value of $z_{\theta}(u)$ as measured in simulations (squares) by analyzing the
statistics of trajectories after the FPT, compared with theoretical values for
different $\theta$. Notice the divergences for small $u$ towards $\pm\infty$,
which enable us to select the value of $\theta$ to minimize these divergences.
(d) Same as (b) for different $H$, with $T=0$ fixed. “Star” symbols are the
simulation results of Ref. Krug1997 .

The second class of problems (type II) corresponds to the determination of
$\theta$ in an idealized situation where a given trajectory for $x(t<0)$ is
assumed to be accessible and observed at all times $t<0$ ; this can be
realized in various settings, ranging from single particle tracking techniques
in the context of transport in complex systems, to the monitoring of the value
of an asset in the context of financial markets. Here we aim at quantifying
the impact of such observation of the system in the past ($t<0$) on its future
dynamics ($t>0$). It is known yaglom1955correlation ; GRIPENBERG1996 ;
anh2004prediction ; inoue2012prediction  that the mean future trajectory (for
$t>0$), conditional to a given observation in the past $x(t<0)$, can be
expressed as a linear combination of all positions in the past. Our approach
makes it possible to determine quantitatively the exponent $\theta$
characterizing the probability $S(t)$ of not crossing this average conditional
trajectory, or of not reaching a fixed threshold above (or below) it, see Fig.
2(a). Strikingly, we find that the value of $\theta$ can be significantly
larger than the value $\theta=1-H$ obtained in absence of any prior
observation of the system. It does not depend on the particular realization of
the observed past trajectory, but only on the fact that this observation is
available. This thus shows that the very observation of the system can
drastically impact the future first-passage statistics, and in fact
effectively accelerate the dynamics at large times because $\theta\geq 1-H$
for all values of $H$, irrespective of the persistent or antipersistent nature
of the process. The results in Fig. 2(b) show again a good agreement between
the predicted values of $\theta$ and simulations. As above a perturbation
expansion of our formalism can be performed for weakly non-Markovian
processes, leading to the explicit result

$\displaystyle\theta_{\mathrm{II}}=1-H+4\ln
2(H-1/2)^{2}+\mathcal{O}((H-1/2)^{3}),$ (9)

which is supported by our simulations [Fig. 2(b)].

Figure 2: Persistence for fBM conditioned on the past trajectory (type II
problem) (a) Definition of the problem: assume that a particular trajectory is
observed for $t<0$ (black curve). The blue line is one realization of
trajectory for $t>0$, the dashed red line represents the average trajectory
given that the past trajectory is observed. The time $t^{*}$ is the first
crossing time to this predicted average trajectory, and the persistence
exponent characterizes the probability that $t^{*}>t$, for large $t$. (b)
Comparison between values of $\theta$ for the fBM conditioned on the past
trajectory, obtained in simulations (symbols), our theoretical approach (black
line) and perturbation expansion (dashed red line). We also indicate the value
$\theta=1-H$ for non-conditioned fBM.

Persistence in higher dimensions.  Our theory can be generalized to the case
of an isotropic Gaussian random process $\mathbf{x}(t)$ evolving in a space of
dimension $d>1$. In this case, to define the survival probability we replace
the condition of reaching a threshold by the condition of reaching a target.
To the best of our knowledge, in this case the persistence exponent has not
been investigated in the literature for non-Markovian walks with non-
stationary initial conditions, despite its obvious relevance to reactivity
problems in complex systems. Here we restrict ourselves to the case where a
target, even point-like, is found with probability one (compact case, when
$dH<1$). It turns out that very few changes are needed to generalize the
theory in $d$ dimensions, generalized versions of the equations are presented
in SI (Section C), and we restrict ourselves to $H>1/(2+d)$. Fig. 3 shows
simulation results when $\mathbf{x}(t)$ is the position of a monomer in
various polymer models: semi-flexible or flexible chains, or fractal
hyperbranched flexible macromolecules. It is found that our theory captures
quantitatively the dependence of the persistence exponents on the temperature
quench for all these models. This dependence on the temperature quench shows
that the exponents describing the kinetics of absorption to a target is
significantly modified by preparing the system with non-stationary initial
conditions. For example, if the random walker is a tagged monomer of a
macromolecule, a non-equilibrium condition could be obtained by a temperature
quench, or by imposing a constraint, such as a geometric confinement or an
external field, that is relaxed at $t=0$. Alternatively, if memory effects of
the random walker come from its interactions with a surrounding viscoelastic
medium, a non-equilibrium initial state could be obtained by imposing a sharp
change of the parameters characterizing this medium. In these cases, we
predict that the reaction kinetics to a target can be deeply impacted, and
display non-trivial exponents quantified by our approach (as soon as the
covariance $\sigma$ can be calculated).

Conclusion. In Fig. 6, simulation data for $\theta$ for all the models
considered in this work are recapitulated and compared with the values
predicted by our approach. The data collapse shows an excellent agreement and
validates our method. The slight departures of simulations from theory occur
only when $2H-\theta>0$ (see Fig. 6), in agreement with our previous remark
that our selection criterium may not be valid anymore in this regime.
Similarly, the curves which are the least precise on Fig. 3 are those for
which $H$ is very close to $1/(2+d)$, where the theory is not expected to give
accurate results anymore. Altogether, this shows that our theory provides a
non perturbative, constructive, quantitative determination of the persistence
exponents for general Gaussian stochastic processes with non-stationary
initial conditions, which typically model the relaxation after an initial
perturbation of systems with non Makovian dynamics, such as tracer particles
in complex environments with many interacting degrees of freedom. Our results
demonstrate that initial perturbations can have a deep, long lived impact on
the first-passage statistics of non-Markovian processes. Importantly, our
theory also predicts non-trivial exponents in dimension higher than one, and
thus opens the way to the quantification and control of reaction kinetics for
complex systems with non-equilibrium initial conditions.

Figure 3: Persistence exponents in dimensions higher than one (a) $d=2$ and
(b) $d=3$. These results are obtained when the random walker is a tagged
monomer in different polymer models: semi-flexible chains ($H=3/8$), flexible
chain without hydrodynamic interactions ($H=1/4$), macromolecule of fractal
architecture (Vicsek fractal, $H\simeq 0.203$). Figure 4: Comparison of
theoretical values of the persistence exponent versus simulations for all data
shown in Figs. 1, 2, 3.

###### Acknowledgements.

Computer time for this study was provided by the computing facilities MCIA
(Mesocentre de Calcul Intensif Aquitain) of the Université de Bordeaux and of
the Université de Pau et des Pays de l’Adour. T-V. M. and T. G. acknowledge
the support of the grant ComplexEncounters, ANR-21-CE30-0020-01.

Supplementary Information

In this Supplementary information, we provide:

  * •

a characterization of Gaussian processes that have stationary increments at
long times (Section A),

  * •

the reason why persistence exponents for such processes cannot be obtained by
standard methods (Section B),

  * •

the list of the stochastic processes (and associated simulation algorithms)
investigated in this work (Section C),

  * •

a detailed description of the derivation of Eqs. (3), (4), (6) of the main
text (Section D), including the description of how $\theta$ is numerically
obtained (Section D.5).

  * •

a discussion on the restrictions of the theory (Section E),

  * •

a perturbative analysis of our theory around the Markovian Brownian motion
($H\to 1/2$, Section F),

  * •

the analysis of the quenched fBm in the limit of infinite initial temperature
(Section G).

##  Appendix A Characterization of Gaussian processes with stationary
increments at long times

We consider a one-dimensional Gaussian random walk $x(t)$ evolving with
continuous time $t\geq 0$. Since the process is Gaussian, it is entirely
defined by its mean $\langle x\rangle$, which we assume to be constant with
time (symmetric walk), and by its covariance
$\mathrm{Cov}(x(t),x(t^{\prime}))=\sigma_{0}(t,t^{\prime})$. At long times, we
assume that $\sigma_{0}(t,t^{\prime})$ admits the scaling behavior:

$\displaystyle\sigma_{0}(t,t^{\prime})\underset{t,t^{\prime}\to\infty}{\sim}t^{2H}\
G(t/t^{\prime})\equiv\sigma(t,t^{\prime}),$ (10)

where $G$ is a scaling function, $H$ is the Hurst exponent, and the symbol
$\sim$ represents the mathematical asymptotic equivalence. The above scaling
behavior is standard (see e.g. Ref. ReviewBray ). Note also that
$\sigma(t,t^{\prime})$ is the covariance function at long time scales and is
defined by Eq. (10). We will restrict our-selves to $H>0$ (so that the
particle is not trapped at long times at any position of space), and also to
$H<1$ (the long-time process is then non-smooth ReviewBray ).

In this work, we study Gaussian stochastic processes that have stationary
increments at long times, thus displaying transient aging. In other words, we
assume that the statistics of increments $x(t+\tau)-x(t)$ do not depend on the
time elapsed $t$, if $t$ is large enough. For such processes, one can define a
“stationary covariance” $\sigma_{s}$ by:

$\displaystyle\sigma_{s}(\tau,\tau^{\prime})=$
$\displaystyle\lim_{t\to\infty}\langle[x(t+\tau)-x(t)][x(t+\tau^{\prime})-x(t)]\rangle.$
(11)

For Gaussian processes whose covariance satisfies Eq. (10), the existence of
$\sigma_{s}$ requires conditions on the behavior of the scaling function $G$
which we determine now. Using Eq. (10), we obtain

$\displaystyle\sigma_{s}(\tau,\tau^{\prime})=$
$\displaystyle\lim_{t\to\infty}\left\\{(t+\tau)^{2H}G\left(\frac{t+\tau^{\prime}}{t+\tau}\right)-t^{2H}\left[G\left(\frac{t}{t+\tau}\right)+G\left(\frac{t}{t+\tau^{\prime}}\right)-G\left(1\right)\right]\right\\}$
(12)

Thus, the quantity $\sigma_{s}$ exists when the behavior of $G(u)$ near $u=1$
is

$\displaystyle G(u)\underset{u\to
1,u<1}{=}G(1)[1+H(1-u)]-\kappa_{s}(1-u)^{2H}+\text{h.o.t},$ (13)

which holds for $u<1$. Here, h.o.t means higher order terms, and $\kappa_{s}$
is a positive constant. For $u>1$ the behavior of $G$ is deduced from the
condition that $\sigma(t,t^{\prime})=\sigma(t^{\prime},t)$ is symmetric. Note
that specifying the linear term (proportional to $1-u$) is not necessary for
$H<1/2$ but it dominates for superdiffusive processes so that it is important
in this case. In fact, all correlators we will study in this work satisfy the
above condition, irrespectively on the value of $H$. With the above condition,
the stationary covariance reads

$\displaystyle\sigma_{s}(\tau,\tau^{\prime})=\kappa_{s}\left(\
|\tau|^{2H}+|\tau^{\prime}|^{2H}-|\tau^{\prime}-\tau|^{2H}\right)=\tau^{2H}G_{s}(\tau/\tau^{\prime}).$
(14)

This means that the increments in the future of $t$ have the same covariance
function as that of the fractional Brownian motion (fBm), multiplied by a
generalized transport coefficient $\kappa_{s}$. Here, $G_{s}$ is the scaling
function associated to the stationary state and is defined by the above
equation. The above formula holds also when $\tau,\tau^{\prime}$ are negative.

Alternatively, to define the stationary property at long times, we could
require that, if we look only at trajectories satisfying the constraint
$x(t)=x^{*}$, and that one looks at the trajectories in the future of this
time $t$, then the statistics of these conditioned trajectories do not depend
on $t$, if $t$ is large enough. In other words, we could define
$\sigma_{s}(\tau,\tau^{\prime})$ as the limiting value of the covariance of
$x(t+\tau)$ and $x(t+\tau^{\prime})$ given that $x(t)=x^{*}$ is fixed, when
$t$ is large enough:

$\displaystyle\sigma_{s}(\tau,\tau^{\prime})=\lim_{t\to\infty}\mathrm{Cov}(x(t+\tau),x(t+\tau^{\prime})|x(t)=x^{*})=\lim_{t\to\infty}\left\\{\sigma_{0}(t+\tau,t+\tau^{\prime})-\frac{\sigma_{0}(t+\tau,t)\sigma_{0}(t,t+\tau^{\prime})}{\sigma_{0}(t,t)}\right\\},$
(15)

where the notation $\text{Cov}(A,B|E)$ represents the covariance of $A,B$
given that the event $E$ is realized, and the second equality is deduced from
general formulas on conditional Gaussian processes (see e.g. Ref. Eaton1983 ).
It turns out that the condition for the existence of $\sigma_{s}$ defined in
(14) or (15) is the same and is given by Eq. (13), and that the two
definitions of $\sigma_{s}$ are identical. Hence, requiring that the
increments of the trajectories are stationary at long time, or imposing a
condition on the stationarity of conditional covariances, is actually the
same. It is instructive to consider the convergence of the conditional
covariance towards the stationary covariance:

$\displaystyle\mathrm{Cov}(x(t+t_{1}),x(t+t_{2})|x(t)=0)-\sigma_{s}(t_{1},t_{2})\underset{t\to\infty}{\sim}-\frac{\kappa_{s}^{2}(t_{1}t_{2})^{2H}}{t^{2H}G(1)}.$
(16)

This convergence is therefore algebraic, which was expected due to the absence
of characteristic relaxation time.

In summary, in this Section we have characterized the Gaussian processes that
have long time stationary increments: these are the Gaussian processes that
display a scaling function $G$ satisfying the condition (13).

##  Appendix B The reason why standard methods to calculate persistence
exponents fail when applied to Gaussian processes with stationary increments
at long times

For one-dimensional Gaussian stochastic processes, a standard method to
calculate the persistence exponent $\theta$ (defined so that $S(t)\propto
1/t^{\theta}$ at long times) consists in applying the Lamperti transform to
obtain a stationary Gaussian process, and then to apply the independent
interval approximation (IIA). This method is not exact but gives in general
very good estimates of persistence exponents Krug1997 . The Lamperti transform
consists in defining logarithmic times $T\equiv\ln t$ and rescaling $x(t)$ by
its root mean square to obtain a stationary process, i.e. defining the
stochastic process $X(T)\equiv x(t)/[G(1)t^{2H}]^{1/2}$. The covariance
function of the transformed process $X(T)$ (still Gaussian) reads

$\displaystyle\langle
X(T_{1})X(T_{2})\rangle=\frac{t_{1}^{2H}G(t_{1}/t_{2})}{G(1)t_{1}^{H}t_{2}^{H}}=e^{-H(T_{2}-T_{1})}\frac{G(e^{-(T_{2}-T_{1})})}{G(1)}\equiv
a(|T_{2}-T_{1}|).$ (17)

This correlation function depends only on the difference $T_{2}-T_{1}$,
meaning that the process $X(T)$ is stationary. For small positive $T$, with
the hypothesis (13), this correlation function behaves as

$\displaystyle a(T)\underset{T\to
0}{\simeq}(1-HT+...)[G(1)(1+HT)-\kappa_{s}T^{2H}+...]/G(1)=1-[\kappa_{s}/G(1)]T^{2H}+\text{h.o.t}.$
(18)

For $0<H<1$, the above expansion means that $X(T)$ is a non-smooth Gaussian
process, because smooth Gaussian processes admit a correlation function
behaving as $a(T)=a(0)-bT^{2}+...$ for $T\to 0$. Now, the IIA method consists
in evaluating the statistics of the time between two successive zeros of a
smooth stationary process, but by construction for a non-smooth process this
very concept of interval between two zeros is ill-defined. Hence by
construction the IIA method can be applied only to smooth Gaussian processes.
Therefore it cannot be applied to evaluate persistence exponents in our
problems, where the dynamics displays stationary increments at long times.
This has been noted in Ref. Krug1997  for the study of persistence in
fluctuating interfaces, but what we have shown here is that this restriction
is more general and associated to the fact that the process has stationary
increments at long times.

##  Appendix C List of stochastic processes investigated in this work

###  C.1 “Quenched” fBm (type I process)

We describe here the stochastic processes for which we will study persistence
exponents. We will consider a process which we call “quenched” fBm, for which
the correlator at long times is

$\displaystyle\sigma(t,t^{\prime})=\sigma_{\text{I}}(t,t^{\prime})\equiv\frac{(1-T)(t+t^{\prime})^{2H}+T[t^{2H}+t^{\prime
2H}]-|t-t^{\prime}|^{2H}}{[(1-T)2^{2H}+2T]}.$ (19)

The function $G$ associated to this correlator satisfies the condition (13),
meaning that it corresponds to a process $x(t)$ that becomes with stationary
increments at long times. Depending on $H$, this correlator corresponds to
various physical models (which we used to simulate $x(t)$ for various values
of $H$). For example, in $d=1$, $x(t)$ can be seen as the local height of a
one dimensional interface $x(t)=h(s=0,t)$ which obeys the dynamics

$\displaystyle\partial_{t}h=-(-\partial_{s}^{2})^{1/(2-4H)}h(s,t)+\xi(s,t),\hskip
28.45274pt\langle\xi(s,t)\xi(s^{\prime},t^{\prime})\rangle=2T_{1}\delta(s-s^{\prime})\delta(t-t^{\prime}),$
(20)

where $T_{1}$ is the temperature during the dynamics. Here, we will focus on
the cases $H=1/4$ (Edwards-Wilkinson dynamics) and $H=3/8$ (Mullins-Herring
dynamics). The correlator (19) is obtained when $x(t)=h(s=0,t)$ is calculated
with an initial condition corresponding to a stationary state of Eq. (20) with
initial temperature $T_{1}=T$, while the dynamics takes place at temperature
$T_{1}=1$. The fact that (19) is obtained after a sharp change of temperature
in a model of collective dynamics justifies the name “quenched fBm”.
Interestingly, the correlator (19) with $T=0$ and $H=1/4$ is also that of a
tagged particle in single file diffusion in crowded narrow channels, for which
the position of all other particles is initially fixed (which corresponds to
the zero temperature state in this problem) krapivsky2015dynamical . Next, if
one replaces $h(s,t)$ by $\mathbf{r}(s,t)$ in Eq. (20), one obtains the
equation for the dynamics of the positions of monomers $\mathbf{r}(s,t)$ at
time $t$ and position $s$ along the chain, for the Rouse chain ($H=1/4$, a
bead-spring chain) or the semi-flexible chain ($H=3/8$). Other values of $H$
can be obtained for the dynamics of the position $\mathbf{r}_{i}(t)$ of
monomers in hyperbranched (fractal) macromolecules

$\displaystyle\partial_{t}\mathbf{r}_{i}=-\sum_{j}A_{ij}\mathbf{r}_{j}(t)+\mathbf{f}(s,t),\hskip
28.45274pt\langle
f_{i,\alpha}(s,t)f_{j,\beta}(s^{\prime},t^{\prime})\rangle=2T_{1}\delta_{ij}\delta(t-t^{\prime})\delta_{\alpha,\beta},$
(21)

where $\alpha,\beta$ are spatial coordinates ($x,y,z$ in $d=3$) and $A_{ij}$
is the dynamical matrix, with $A_{ij}=-1$ if beads $i,j$ are connected and
zero in the contrary, and $A_{ii}$ is the number of connected neighbors around
bead $i$.

Simulations.-  In practice, we have simulated processes $x(t)$ (or
$\mathbf{r}(t)$ in $d$ dimensions) by integrating numerically a discretized
versions of Eq. (20) (for $H=1/4$ and $H=3/8$) or Eq. (21), for which we took
the connectivity matrix of a Vicsek fractal of functionality $f=4$
blumen2004generalized , with $H=\ln(3)/\\{2\ln[3(1+f)]\\}\simeq 0.203$.
Initially the polymer chain (or interface) was prepared at temperature $T$
with a reactive monomer at distance $x_{0}$ from the target and we recorded
the time to the crossing of $x=0$ (in $d=1$), or the time to reach a sphere of
radius $a$ in $d$ dimensions. To simulate the dynamic of VFs, we have
decomposed the dynamics on eigenmodes

$\displaystyle\mathbf{r}_{n}=\sum_{i}b_{i}\mathbf{a}_{i}\ ;\
\partial_{t}\mathbf{a}_{i}=-\lambda_{i}\mathbf{a}_{i}+\mathbf{f}_{i}(t),\langle
f_{i}(t)f_{j}(t^{\prime})\rangle=2T_{1}\delta(t-t^{\prime})\delta_{ij}$ (22)

where $\mathbf{r}_{n}$ is the reactive monomer, $\lambda_{i}$ are the
eigenvalues of $A$ and $b_{i}$ the projection of the vector
$(0,...,0,1,0,...0)$ on the eigenspace associated to $\lambda_{i}$. We used
the methods of Ref dolgushev2015contact  to identify iteratively the set of
$(\lambda_{i},b_{i})$. We simulated (22), which is much more efficient than
using (21) since the number of distinct eigenvalues for fractals is much
smaller than the number of beads. In all simulations, we measured $\theta$ by
plotting the survival probability in log-log scale, including the largest
times allowed by the number of accumulated stochastic trajectories. Note that,
in particular for low values of $\theta$, the regime of low values of $S(t)$
could not be reached so that the measurement of $\theta$ is not very accurate.
We checked that the measured values do not depend on microscopic parameters
(initial distance to target, target size, time steps, number of monomers…).

###  C.2 fBm constrained on its past (type II process)

We consider here the statistics of trajectories of a fractional Brownian
motion conditioned on its past value. This process is obtained by considering
a standard fractional Brownian motion, and asking for the statistics of $x(t)$
for $t>0$, given that the past trajectory $x(\tau)$ is known exactly for all
past times $\tau<0$. To understand how this process is built, let us consider
the case that the process is conditioned on $n$ past times $\tau_{i}<0$
(instead of all past times). In this case, general formulas on Gaussian
processes (see e.g. chapter 3 in Ref Eaton1983 ) show that the statistics of
the conditioned process can be obtained from the knowledge of the covariance
matrix $S_{ij}=\sigma_{s}(\tau_{i},\tau_{j})$ (for the unconditioned process,
at the times $\tau_{i}$ at which one imposes the conditions):

$\displaystyle\mathbb{E}(x(t)|[x(\tau_{1}),...x(\tau_{n})])=\sum_{i,j=1}^{n}\sigma_{s}(t,\tau_{i})\
(S^{-1})_{ij}\ x(\tau_{j}),$ (23)
$\displaystyle\mathrm{Cov}(x(t),x(t^{\prime})|[x(\tau_{1}),...x(\tau_{n})])=\sigma_{s}(t,t^{\prime})-\sum_{i,j=1}^{n}\sigma_{s}(t,\tau_{i})\
(S^{-1})_{ij}\ \sigma_{s}(t^{\prime},\tau_{j}).$ (24)

The conditional means and covariances can thus be found at the cost of solving
a linear equation to find the inverse matrix $S^{-1}$. In the continuous limit
$n\to\infty$, the problem of finding $S^{-1}$ is replaced by an integral
equation which was solved by Yaglom (for $H<1/2$) yaglom1955correlation  and
Gripenberg and Norros (for $H>1/2$) GRIPENBERG1996  (see also
anh2004prediction ; inoue2012prediction ), with the result:

$\displaystyle\mathbb{E}(x(t)|[x(\tau),\tau<0])=\begin{cases}\frac{\cos(\pi
H)}{\pi}\int_{0}^{\infty}d\tau\left(\frac{t}{\tau}\right)^{H+1/2}\frac{x(-\tau)}{t+\tau}&(H<1/2),\\\
x(0)+\frac{\cos\pi H}{\pi}\int
dx(-\tau)\int_{0}^{t/\tau}dv\frac{v^{H-1/2}}{v+1}&(H>1/2).\end{cases}$ (25)

Note that in the second line the integration runs over stochastic
trajectories, see Ref. GRIPENBERG1996 . It is clear from Eqs. (23), (24) that
the formula for the conditional covariance is the same as that for the
conditional mean (up to an additive term $\sigma_{s}$) if we replace the past
trajectory $x(\tau_{i})$ by $\sigma_{s}(t^{\prime},\tau_{i})$, so that

$\displaystyle\sigma(t,t^{\prime})=\text{Cov}(x(t),x(t^{\prime})|[x(\tau),\tau<0])=\sigma_{s}(\tau,\tau^{\prime})-\frac{\cos(\pi
H)}{\pi}\times\begin{cases}\int_{0}^{\infty}d\tau\left(\frac{t}{\tau}\right)^{H+1/2}\frac{\sigma_{s}(t^{\prime},-\tau)}{t+\tau}&(H<\frac{1}{2})\\\
\int_{0}^{\infty}d\tau\int_{0}^{t/\tau}dv\frac{v^{H-1/2}}{v+1}\frac{d}{d\tau}\sigma_{s}(t^{\prime},-\tau)&(H>\frac{1}{2})\end{cases}$
(26)

Note that, in the case $H>1/2$ in Eq. (25) we have considered that the
integration over the trajectory can be replaced by a sum over
$x(t_{i+1})-x(t_{i})$ in discrete time, applying the same linear operator on
$\sigma_{s}$ leads to the derivative of $\sigma_{s}$ with respect to $\tau$ in
the above equation. Replacing $\sigma_{s}$ by its value (14), we obtain the
following formula for the covariance of a fBm conditioned on its past:

$\displaystyle\sigma(t,t^{\prime})=\sigma_{\text{II}}(t,t^{\prime})\equiv\kappa_{s}\left\\{-|t-t^{\prime}|+t^{2H}+t^{\prime
2H}+\int_{0}^{\infty}dx\frac{\cos(\pi H)\
[(t+xt^{\prime})^{2H}-t^{2H}-(t^{\prime}x)^{2H}]}{\pi\
x^{H+1/2}(1+x)}\right\\},$ (27)

which is valid both for $H>1/2$ and $H<1/2$. Note that the normalization
$\sigma(t,t)=t^{2H}$ is obtained if we chose $\kappa_{s}=\Gamma(1-H)\
\Gamma(1/2+H)/[4^{H}\sqrt{\pi}]$

Simulations. To simulate a fBm constrained on its past, we used the Hosking
algorithm hosking1984modeling ; Dieker2004  in which the fBm is iteratively
sampled at fixed time steps $t_{i}=n\Delta t$, each position $x(t_{i+1})$
being generated from $x(t_{0}=0),x(t_{1}),...,x(t_{i})$. To take into account
the fact that the past trajectory is known, we defined a time $t_{m}$ and
replaced all positions $x(t_{i})$ with $t_{i}<t_{m}$ by a starting position
$x_{0}$. The next steps $x(t_{j})$ with $t_{j}>t_{m}$ were generated according
to the standard Hosking algorithm. In the limit $t_{m}\to\infty$, this
corresponds to a fBm constrained on its past trajectory after removing its
average value to obtain a symmetric (centered) process.

##  Appendix D Detailed derivation of the equations of the formalism

###  D.1 The link between the covariance of trajectories in the late future of
the first-passage time and the persistence exponent [Derivation of Eq. (3)]

Here, we consider a stochastic process in $d$ dimensions,
$\mathbf{x}(t)=(x(t),y(t),...)$. We suppose it is Gaussian, so that it is
defined by its mean value and its covariance. We also assume that it is
isotropic, so that the mean value $\langle\mathbf{x}(t)\rangle$ is constant
with time, and the covariance between coordinates takes the form
$\text{Cov}[x_{\alpha}(t),x_{\beta}(t^{\prime})]=\delta_{\alpha\beta}\sigma_{0}(t,t^{\prime})$,
where $\alpha,\beta$ label the different spatial coordinates. This is due to
the fact that $\text{Cov}[x_{\alpha}(t),x_{\beta}(t^{\prime})]$ is, at fixed
$t,t^{\prime}$, an isotropic rank-2 tensor and is therefore proportional to
$\delta_{\alpha\beta}$. We further assume that each coordinate $x_{\alpha}(t)$
is a Gaussian stochastic process that has long time stationary increments, so
that $\sigma_{0}$ satisfies the scaling behavior (10) with the condition (13).
If $d=1$, we consider the first-passage problem to a fixed threshold, for
which the persistence $S(t)$, i.e. the probability of not having crossed this
threshold, decays as $S(t)\propto 1/t^{\theta}$. The natural generalization of
this definition to $d\geq 1$ is to consider the first-passage problem to a
target region (which is reduced to a point for $d=1$); here we assume that the
process is compact ($dH<1$), so that the target can be taken as punctual
without any ambiguity. We also assume the process is non-smooth. Our starting
point is a generalized form of the renewal equation

$\displaystyle p(\mathbf{0},t)=\int_{0}^{t}d\tau\
f(\tau)p(\mathbf{0},t|\mathrm{FPT}=\tau),$ (28)

where $p(\mathbf{x},t)$ is the probability density to observe the position
$\mathbf{x}$ (in $d$ dimensions) at time $t$, starting from the initial
conditions at $t=0$, and $f(t)=-\partial_{t}S(t)$ is the probability density
that the FPT is equal to $t$. Here we have used the hypothesis that the
process is rough, and that it is compact, so that a punctual target can be
found by the random walker. The above equation is simply obtained by
partitioning the event of being inside the target at time $t$ over the value
of the first-passage time.

We define $p_{\pi}(x,t)$ as the probability density to observe $\mathbf{x}$ at
a time $t$ after the FPT:

$\displaystyle p_{\pi}(\mathbf{x},t)\equiv\int_{0}^{\infty}d\tau\
p(\mathbf{x},t+\tau|\mathrm{FPT}=\tau)f(\tau).$ (29)

Next, we consider a fixed time $T$, and integrate the renewal equation (28)
between $t=0$ and $t=T$:

$\displaystyle\int_{0}^{T}dt\ p(\mathbf{0},t)=\int_{0}^{T}dt\int_{0}^{t}d\tau
f(\tau)p(\mathbf{0},t|\mathrm{FPT}=\tau)=\int_{0}^{T}d\tau\int_{\tau}^{T}dtf(\tau)p(\mathbf{0},t|\mathrm{FPT}=\tau),$
(30)

where we have simply inverted the order of integration in the second equality.
Setting $t=u+\tau$, we obtain

$\displaystyle\int_{0}^{T}dt\
p(\mathbf{0},t)=\int_{0}^{T}d\tau\int_{0}^{T-\tau}duf(\tau)p(\mathbf{0},\tau+u|\mathrm{FPT}=\tau)=\int_{0}^{T}du\int_{0}^{T-u}d\tau
f(\tau)p(\mathbf{0},\tau+u|\mathrm{FPT}=\tau),$ (31)

where, again the order of integration has been inverted in the second
equality. Now, we integrate Eq. (29) over $t$ between $t=0$ and $t=T$ and use
Eq. (31) to obtain

$\displaystyle\int_{0}^{T}dt\
[p_{\pi}(\mathbf{0},t)-p(\mathbf{0},t)]=\int_{0}^{T}du\int_{T-u}^{\infty}d\tau
f(\tau)p(\mathbf{0},\tau+u|\mathrm{FPT}=\tau).$ (32)

We will evaluate both sides of this equation for large $T$. To evaluate the
right hand side (rhs) term, we assume that

$\displaystyle p(\mathbf{0},\tau+u|\mathrm{FPT}=\tau)\sim
B(\tau/u)/u^{dH}\hskip 28.45274pt(u\to\infty,\tau\to\infty),$ (33)

where $B$ is a scaling function that is assumed to be finite and admit finite
limits for small and large arguments. The above hypothesis is reasonable,
because in the future $u$ of the FPT, one expects that the probability
distribution of $\mathbf{x}$ will be some function whose extension is
proportional $L\sim u^{H}$, so that the probability to be at
$\mathbf{x}=\mathbf{0}$ is $\sim 1/L^{d}=1/u^{dH}$. The prefactor that can
depend on $u/\tau$, and the hypothesis that $B$ is finite is well supported by
our simulations (see e.g. Fig. 5). Now, by definition of $\theta$, we have for
large times

$\displaystyle f(t)\underset{t\to\infty}{\simeq}f_{0}/t^{\theta+1},$ (34)

where $f_{0}$ is a prefactor. The term on the right hand side of Eq. (32) can
be evaluated in the large $T$ limit by setting $\tilde{u}=u/T$ and
$\tilde{\tau}=\tau/T$ and taking the large $T$ limit in the resulting
expression, leading to

$\displaystyle\int_{0}^{T}du\int_{T-u}^{\infty}d\tau
f(\tau)p(\mathbf{0},\tau+u|\mathrm{FPT}=\tau)$
$\displaystyle=T^{2}\int_{0}^{1}d\tilde{u}\int_{1-\tilde{u}}^{\infty}d\tilde{\tau}f(T\tilde{\tau})p(\mathbf{0},T(\tilde{\tau}+\tilde{u})|\mathrm{FPT}=T\tilde{\tau})$
(35) $\displaystyle\underset{T\to\infty}{\simeq}T^{1-\theta-
dH}\int_{0}^{1}d\tilde{u}\ Q(\tilde{u}),$ (36)

where $Q$ is given by

$\displaystyle Q(\tilde{u})$
$\displaystyle=\int_{1-\tilde{u}}^{\infty}d\tilde{\tau}\frac{f_{0}B(\tilde{\tau}/\tilde{u})}{\tilde{\tau}^{1+\theta}\tilde{u}^{dH}}=\frac{f_{0}}{\tilde{u}^{\theta+dH}}\int_{(1-\tilde{u})/\tilde{u}}^{\infty}d\xi\frac{B(\xi)}{\xi^{1+\theta}}.$
(37)

The existence of $Q$ is guaranteed by the fact that $\theta>0$. Its
asymptotics read

$\displaystyle Q(\tilde{u}\to
0)\simeq\frac{f_{0}}{\tilde{u}^{\theta+dH}}\frac{B(\infty)}{\theta\tilde{u}^{-\theta}}\simeq\frac{f_{0}B(\infty)}{\theta\
\tilde{u}^{dH}},\hskip 56.9055ptQ(\tilde{u}\to 1)\simeq
f_{0}\frac{B(0)}{\theta\ (1-\tilde{u})^{\theta}}.$ (38)

The above relations suffice to establish that $\int_{0}^{1}dxQ(x)$ exists if
$\theta<1$ and $dH<1$. The condition $dH<1$ was already assumed (compact
process), and from now on we also assume $\theta<1$ (meaning that the mean
first-passage time is infinite). In these conditions, we deduce the scaling

$\displaystyle\int_{0}^{T}du\int_{T-u}^{\infty}d\tau
f(\tau)p(\mathbf{0},\tau+u|\mathrm{FPT}=\tau)\propto T^{1-\theta-dH}.$ (39)

We now evaluate the term in the left hand side of Eq. (32). At leading order,
dimensional analysis indicate that

$\displaystyle
p(\mathbf{0},t)\underset{t\to\infty}{\sim}\frac{K_{0}}{t^{dH}},\hskip
28.45274ptp(\mathbf{0},t)\underset{t\to\infty}{\sim}\frac{K_{\pi}}{t^{dH}},$
(40)

where $K_{0}$ and $K_{\pi}$ are constants. Since $dH<1$, if $K_{0}\neq
K_{\pi}$, we have at leading order

$\displaystyle\int_{0}^{T}dt\
[p_{\pi}(\mathbf{0},t)-p(\mathbf{0},t)]\underset{T\to\infty}{\sim}\frac{(K_{\pi}-K_{0})}{1-dH}T^{1-dH}.$
(41)

This behavior must be matched with that of the rhs of (32), which behaves as
$T^{1-dH-\theta}$ from Eq. (39). Since $\theta>0$, we obtain a contradiction,
which means that at long times, $p(\mathbf{0},t)$ and $p_{\pi}(\mathbf{0},t)$
are equal at leading order:

$\displaystyle K_{0}=K_{\pi}.$ (42)

We thus postulate an exponent $\alpha>dH$ for the next-to-leading order
correction, such that

$\displaystyle
p_{\pi}(\mathbf{0},t)-p(\mathbf{0},t)\simeq\frac{K_{1}}{t^{\alpha}}+...\hskip
28.45274pt(t\to\infty).$ (43)

In order to find the value of $\alpha$, we have to discuss according to the
value of $\theta$. If $\theta<1-dH$, we see from Eq. (39) that the rhs of Eq.
(32) diverges with $T$ for large $T$. This means that $\alpha<1$, so that the
left hand side term of Eq. (32) also diverges with $T$, as $T^{1-\alpha}$.
Equating the divergences leads to

$\displaystyle\alpha=dH+\theta.$ (44)

In the opposite case $\theta>1-dH$, we see from Eq. (39) that the rhs of Eq.
(32) tends to $0$ with $T$ for large $T$. This means that the following
condition must hold exactly:

$\displaystyle\int_{0}^{\infty}dt\
[p_{\pi}(\mathbf{0},t)-p(\mathbf{0},t)]=0\hskip 28.45274pt(\text{if
}1>\theta>1-dH).$ (45)

Then, using the trick
$\int_{0}^{T}(...)=\int_{0}^{\infty}(...)-\int_{T}^{\infty}(...)$, we find

$\displaystyle-\int_{T}^{\infty}[p_{\pi}(\mathbf{0},t)-p(\mathbf{0},t)]\propto
T^{1-dH-\theta},$ (46)

which leads to $\alpha=dH+\theta$ again. To proceed further, we define the
covariance in the future of the first-passage time,

$\displaystyle\sigma_{\pi}(t,t^{\prime})=\text{Cov}[x(\mathrm{FPT}+t),x(\mathrm{FPT}+t^{\prime})],$
(47)

and we define an exponent $\beta$ characterizing the divergence of the mean
square displacement (MSD) after the FPT as

$\displaystyle\sigma_{\pi}(t,t)-\sigma_{0}(t,t)\simeq c_{1}\ t^{\beta}\hskip
28.45274pt(t\rightarrow\infty).$ (48)

If we assume that the stochastic process in the future of the FPT is Gaussian
(which is well supported by our simulations, see e.g. Fig. 5), we obtain
$p_{\pi}(\mathbf{0},t)\propto\sigma_{\pi}(t,t)^{-d/2}$. This leads to
$\alpha=(d+2)H-\beta$, and thus

$\displaystyle\beta=2H-\theta$ (49)

Inserting Eq. (49) into Eq. (48), we obtain Eq. (3) in the main text. In
summary, here we have identified a link between the persistence exponent
$\theta$ and the exponent $\beta$ which characterizes the MSD of the
trajectories in the future of the FPT at long times. Let us remind our
hypotheses: (i) $dH<1$ (compact process), (ii) $\theta<1$ (infinite mean FPT),
(iii) that the scaling function $B$ characterizing the propagator in the
future of the FPT is finite for small and large arguments, and (iv) the
trajectories in the future of the FPT display Gaussian statistics. The
validity of (iii) and (iv) can be appreciated on Fig. 5 in one example.

Figure 5: Results obtained for type I process, for $H=3/8$, $T=10$ and $d=1$.
For these parameters, we have simulated numerically Eq. (20), for which
$x(t)=h(s=0,t)$ starts at $x_{0}=1$. We have recorded the statistics of the
trajectories after the first-passage to $x=0$, and identified the conditional
variance of $x(t+\mathrm{FPT})$ conditioned on the value of the
$\mathrm{FPT}=\tau$. This quantity is represented in (a) in rescaled
variables; if the Gaussian hypothesis holds then the function appearing $B$ in
Eq. (33) is $B(y)=1/\sqrt{2\pi\sigma_{\pi}(1,1|y)}$, hence these data support
the fact that $B(y)$ is finite for $y\to\infty$ and $y\to 0$. The black line
represents the ansatz $B$ in Eq. (82). In (b) we have represented the
normalized distribution of $x(t+\text{FPT})$ for various values of $t$, small
deviations from Gaussian behavior can be observed at small $t$ and disappear
at larger values.

###  D.2 Self-consistent equations for the covariance of trajectories in the
future of the FPT

We now write the generalization of the renewal equation (28) for the
probability to observe $\mathbf{x}=0$ at $t$ and the sequence of positions
$\mathbf{x}_{1},\mathbf{x}_{2}$ at times $t+t_{1},t+t_{2}$, where
$t_{1},t_{2}$ are fixed times:

$\displaystyle
p(\mathbf{0},t;\mathbf{x}_{1},t+t_{1};\mathbf{x}_{2},t+t_{2})=\int_{0}^{t}d\tau
f(\tau)\
p(\mathbf{0},t;\mathbf{x}_{1},t+t_{1};\mathbf{x}_{2},t+t_{2}|\mathrm{FPT}=\tau).$
(50)

We also define the probability density to observe the sequence
$(\mathbf{0};\mathbf{x}_{1};\mathbf{x}_{2})$ at respective times
$t,t+t_{1},t+t_{2}$ in the future of the FPT, averaged over the FPT
distribution:

$\displaystyle
p_{\pi}(\mathbf{0},t;\mathbf{x}_{1},t+t_{1};\mathbf{x}_{2},t+t_{2})=\int_{0}^{\infty}d\tau
f(\tau)p(\mathbf{0},t+\tau;\mathbf{x}_{1},t+\tau+t_{1};\mathbf{x}_{2},t+\tau+t_{2}|\mathrm{FPT}=\tau).$
(51)

Repeating the steps that had led to Eq. (32), we obtain

$\displaystyle\int_{0}^{T}dt[p_{\pi}(\mathbf{0},t;\mathbf{x}_{1},t+t_{1};\mathbf{x}_{2},t+t_{2})$
$\displaystyle-p(\mathbf{0},t;\mathbf{x}_{1},t+t_{1};\mathbf{x}_{2},t+t_{2})]=$
$\displaystyle\int_{0}^{T}du\int_{T-u}^{\infty}d\tau
f(\tau)p(\mathbf{0},\tau+u;\mathbf{x}_{1},u+\tau+t_{1};\mathbf{x}_{2},u+\tau+t_{2}|\mathrm{FPT}=\tau).$
(52)

Here we assume, without loss of generality, that the process is centered (zero
average) ; we stress that this has no consequence on the value of $\theta$
which does not depend on microscopic details. We multiply the above equation
by $x_{1}x_{2}$ and integrate over all positions to obtain the equation for
the covariances:

$\displaystyle\int_{0}^{T}dt\ [$ $\displaystyle
p_{\pi}(\mathbf{0},t)\sigma_{\pi}^{*}(t_{1},t_{2}|t)-p(\mathbf{0},t)\sigma^{*}(t_{1},t_{2}|t)]=\int_{0}^{T}du\int_{T-u}^{\infty}d\tau
f(\tau)p(\mathbf{0},\tau+u|\mathrm{FPT}=\tau)\Sigma(t_{1},t_{2},u,\tau)\equiv
C(t_{1},t_{2},T).$ (53)

Here, the term $C$ is defined by this equation, and we have used the following
notations for conditional covariances:

$\displaystyle\sigma^{*}(t_{1},t_{2}|t)=\mathrm{Cov}(x(t+t_{1}),x(t+t_{2})|x(t)=0),$
(54)
$\displaystyle\sigma_{\pi}^{*}(t_{1},t_{2}|t)=\mathrm{Cov}(x(\mathrm{FPT}+t+t_{1}),x(\mathrm{FPT}+t+t_{2})|x(\mathrm{FPT}+t)=0),$
(55)
$\displaystyle\Sigma(t_{1},t_{2},u,\tau)=\text{Cov}(x(u+\tau+t_{1}),x(u+\tau+t_{2})|x(u+\tau)=0,\text{FPT}=\tau).$
(56)

Now, the main problem consists in evaluating $\Sigma$, which is difficult
since it depends it-self on the value of the first-passage time. Nevertheless,
we argue below that it is reasonable to assume that, when a long time has
ocured after the first-passage, this conditional covariance becomes equal to
the long-time stationary covariance:

$\displaystyle\Sigma(t_{1},t_{2},u,\tau)\sim\sigma_{s}(t_{1},t_{2})\hskip
56.9055pt(\tau,u\to\infty).$ (57)

Here, we consider the above equality as an hypothesis, whose validity will be
discussed later (see Section E). If Eq. (57) holds, then the term $C$ can be
evaluated in the limit $T\to\infty$ as

$\displaystyle
C(t_{1},t_{2},T)=T^{2}\int_{0}^{1}d\tilde{u}\int_{1-\tilde{u}}^{\infty}d\tilde{\tau}f(T\tilde{\tau})p(\mathbf{0},T[\tilde{\tau}+\tilde{u}|\mathrm{FPT}=T\tilde{\tau})\Sigma(t_{1},t_{2},\tilde{u}T,\tilde{\tau}T)$
$\displaystyle\underset{T\to\infty}{\sim}T^{1-\theta-
dH}\int_{0}^{1}d\tilde{u}\int_{1-\tilde{u}}^{\infty}d\tilde{\tau}\frac{f_{0}B(\tilde{\tau}/\tilde{u})}{\tilde{\tau}^{1+\theta}\tilde{u}^{dH}}\sigma_{s}(t_{1},t_{2}),$
(58)

with the same notations as in Eq. (36). What is important here is the scaling
with $T$. We see that if $\theta>1-dH$, $C$ tends to zero in the large $T$
limit, this suggests that in the large $T$ limit, Eq. (53) becomes

$\displaystyle\int_{0}^{\infty}dt\
[p_{\pi}(\mathbf{0},t)\sigma_{\pi}^{*}(t_{1},t_{2}|t)-p(\mathbf{0},t)\sigma^{*}(t_{1},t_{2}|t)]=0.$
(59)

This equation can be written explicitly as a function of the covariances
$\sigma_{0}$ and $\sigma_{\pi}$:

$\displaystyle\mathcal{H}(\tau,\tau^{\prime})\equiv\int_{0}^{\infty}dt\Bigg{\\{}\frac{\left[\sigma_{\pi}(t+\tau,t+\tau^{\prime})-\frac{\sigma_{\pi}(t+\tau,t)\sigma_{\pi}(t+\tau^{\prime},t)}{\sigma_{\pi}(t,t)}\right]}{\sigma_{\pi}(t,t)^{d/2}}-\frac{\left[\sigma_{0}(t+\tau,t+\tau^{\prime})-\frac{\sigma_{0}(t+\tau,t)\sigma_{0}(t+\tau^{\prime},t)}{\sigma_{0}(t,t)}\right]}{\sigma_{0}(t,t)^{d/2}}\Bigg{\\}}=0.$
(60)

This equation defines the covariance matrix for $\theta>1-dH$. However, it
cannot be true for $\theta<1-dH$ since the above integral does not exist. Let
us consider now the case $\theta<1-dH$, for which we see that $C$ diverges for
large $T$ [Eq. (58)]. Our approach consists in combining Eqs. (53) and Eq (32)
[multiplied by $\sigma_{s}(t_{1},t_{2})$] to obtain

$\displaystyle\int_{0}^{T}dt\\{p_{\pi}(\mathbf{0},t)[\sigma_{\pi}^{*}(t_{1},t_{2},t)-$
$\displaystyle\sigma_{s}(t_{1},t_{2})]-p(\mathbf{0},t)[\sigma^{*}(t_{1},t_{2}|t)-\sigma_{s}(t_{1},t_{2})]\\}$
$\displaystyle=$ $\displaystyle\int_{0}^{T}du\int_{T-u}^{\infty}d\tau
f(\tau)p(\mathbf{0},\tau+u|\mathrm{FPT}=\tau)[\Sigma(t_{1},t_{2},u,\tau)-\sigma_{s}(t_{1},t_{2})]\equiv
R(t_{1},t_{2},T),$ (61)

where $R$ is defined as the rhs of this equation. We note that $R$ is no
longer of order $T^{1-dH-\theta}$, this problematic term has been removed.
Assuming that $R$ vanishes for large $T$ (this will be discussed in Section
E), we obtain:

$\displaystyle\mathcal{H}(\tau,\tau^{\prime})\equiv\int_{0}^{\infty}dt\Bigg{\\{}\frac{1}{\sigma_{\pi}(t,t)^{d/2}}\Bigg{[}$
$\displaystyle\sigma_{\pi}(t+\tau,t+\tau^{\prime})-\frac{\sigma_{\pi}(t+\tau,t)\sigma_{\pi}(t+\tau^{\prime},t)}{\sigma_{\pi}(t,t)}-\sigma_{s}(\tau,\tau^{\prime})\Bigg{]}$
$\displaystyle-\frac{1}{\sigma_{0}(t,t)^{d/2}}\Bigg{[}$
$\displaystyle\sigma_{0}(t+\tau,t+\tau^{\prime})-\frac{\sigma_{0}(t+\tau,t)\sigma_{0}(t+\tau^{\prime},t)}{\sigma_{0}(t,t)}-\sigma_{s}(\tau,\tau^{\prime})\Bigg{]}\Bigg{\\}}=0.$
(62)

This equation is valid for $\theta<1-dH$; in fact it is also correct for
$\theta>1-dH$, as can be seen by combining Eqs. (60) and (45). For
$\theta=1-dH$ and $d=1$, we recover the equation for processes with stationary
increments guerin2016mean . In summary, here we have obtained self-consistent
equations [Eqs. (60) and (62)] for the covariance in the future of the first-
passage time for non-stationary Gaussian processes.

###  D.3 Equation for the large time behavior of the trajectories after the
first-passage [Derivation of Eq. (4)]

Up to now, the equations for $\sigma_{\pi}(t,t^{\prime})$ involve all time
scales of the process, including microscopic ones. Here, our goal is to
isolate the contribution of large times. First, we note that it is natural to
assume that $\sigma_{\pi}(t,t^{\prime})\sim t^{2H}G_{\pi}(t/t^{\prime})$ for
large $t,t^{\prime}$. We see that $\mathcal{H}(\tau,\tau^{\prime})$ in Eqs.
(60) and (62) is proportional to $\tau^{1-(d-2)H}$ multiplied by a scaling
function of $\tau/\tau^{\prime}$, and the fact that this term must vanish
imposes that $G_{\pi}(x)=G(x)$ for all $x$. This equality generalizes the
result (42). Therefore, we look at large time corrections for $\sigma_{\pi}$,
of the type

$\displaystyle\sigma_{\pi}(t,t^{\prime})\simeq\sigma_{0}(t,t^{\prime})+\rho(t,t^{\prime})+...\hskip
28.45274pt\rho(t,t^{\prime})=t^{2H-\theta}z(t/t^{\prime}),$ (63)

with $z$ a scaling function. In the following, we look for an equation that
defines the function $z(x)$.

First, we consider the case $\theta<1-dH$, for which $\sigma_{\pi}$ satisfies
Eq. (62), which we write for $\tau=Tv$ and $\tau^{\prime}=Tv^{\prime}$:

$\displaystyle\mathcal{H}(Tv,Tv^{\prime})=\int_{0}^{\infty}du\Bigg{\\{}$
$\displaystyle\frac{T}{\sigma_{\pi}(Tu,Tu)^{d/2}}\left[\sigma_{\pi}(T(u+v),T(u+v^{\prime}))-\frac{\sigma_{\pi}(T(u+v),Tu)\sigma_{\pi}(T(u+v^{\prime}),Tu)}{\sigma_{\pi}(Tu,Tu)}-\sigma_{s}(Tv,Tv^{\prime})\right]$
$\displaystyle-\frac{T}{\sigma_{0}(Tu,Tu)^{d/2}}$
$\displaystyle\left[\sigma_{0}(T(u+v),T(u+v^{\prime}))-\frac{\sigma_{0}(T(u+v),Tu)\sigma_{0}(T(u+v^{\prime}),Tu)}{\sigma_{0}(Tu,Tu)}-\sigma_{s}(Tv,Tv^{\prime})\right]\Bigg{\\}}=0,$
(64)

where we have set $t=Tu$. Using (63), the evaluation of all terms when
$T\to\infty$ leads to

$\displaystyle\mathcal{H}$
$\displaystyle(Tv,Tv^{\prime})\underset{T\to\infty}{\sim}T^{1-dH+2H-\theta}\int_{0}^{\infty}\frac{du}{\sigma(u,u)^{d/2}}\Bigg{\\{}\rho(u+v,u+v^{\prime})-\frac{\sigma(u+v,u)\rho(u+v^{\prime},u)}{\sigma(u,u)}-\frac{\sigma(u+v^{\prime},u)\rho(u+v,u)}{\sigma(u,u)}$
$\displaystyle+\rho(u,u)\frac{\sigma(u+v,u)\sigma(u+v^{\prime},u)}{\sigma(u,u)^{2}}-\frac{d\
\rho(u,u)}{2\
\sigma(u,u)}\left[\sigma(u+v,u+v^{\prime})-\frac{\sigma(u+v,u)\sigma(u+v^{\prime},u)}{\sigma(u,u)}-\sigma_{s}(v,v^{\prime})\right]\Bigg{\\}}=0.$
(65)

It can be checked that this integral exists for $\theta<1-dH$. In particular,
the integrand of the above expression for small $u$ is
$\sim[\sigma(v,v^{\prime})-\sigma_{s}(v,v^{\prime})]/u^{dH+\theta}$, which is
integrable since we have $\theta<1-dH$. Since $\mathcal{H}$ must vanish at all
orders, equating the above expression to zero gives Eq. (4) of the main text,
that defines $\rho(v,v^{\prime})$ for $\theta<1-dH$.

If $\theta>1-dH$, we follow the same steps as for $\theta<1-dH$, but we start
from Eq. (60) [instead of (62)]. We obtain

$\displaystyle\mathcal{H}$
$\displaystyle(Tv,Tv^{\prime})\underset{T\to\infty}{\sim}T^{1-dH+2H-\theta}\int_{0}^{\infty}\frac{du}{\sigma(u,u)^{d/2}}\Bigg{\\{}\rho(u+v,u+v^{\prime})-\frac{\sigma(u+v,u)\rho(u+v^{\prime},u)}{\sigma(u,u)}-\frac{\sigma(u+v^{\prime},u)\rho(u+v,u)}{\sigma(u,u)}$
$\displaystyle+\rho(u,u)\frac{\sigma(u+v,u)\sigma(u+v^{\prime},u)}{\sigma(u,u)^{2}}-\frac{d\
\rho(u,u)}{2\
\sigma(u,u)}\left[\sigma(u+v,u+v^{\prime})-\frac{\sigma(u+v,u)\sigma(u+v^{\prime},u)}{\sigma(u,u)}-\sigma(v,v^{\prime})\right]\Bigg{\\}}=0.$
(66)

Note that the only difference with Eq. (65) is the last term, where $\sigma$
appears in place of $\sigma_{s}$. For $u\to\infty$, the integrand of the above
expression is $\sim
d[\sigma_{s}(v,v^{\prime})-\sigma(v,v^{\prime})]/(2u^{dH+\theta})$, which is
integrable since $\theta>1-dH$. Equating the above expression to zero leads to
Eq. (4) in the main text for $\theta>1-dH$.

###  D.4 The equation for $z(x)$ [Derivation of Eq. (6)]

Due to scaling properties, the solutions of Eqs. (65) and (66) can be found
under the scaling form $\rho(v,v^{\prime})=v^{2H-\theta}z(v/v^{\prime})$,
where the single-valued function $z(x)$ satisfies a linear equation which we
determine now. Using the same arguments as in Section A, we can show that, if
the process in the future of the FPT has long-time stationary increments, and
if $z(x)$ has a linear term in its expansion near $x=1$, then

$\displaystyle z(x)\underset{x\to 1}{\simeq}z(1)[1+\beta/2(1-x)]+...$ (67)

Since it is reasonable to believe that the process in the future of the FPT
admits stationary increments at long times, it is also reasonable to assume
the above equality holds, and it is actually the case for the functions $z$
determined at perturbative order, see Eqs. (108),(129). This suggests to pose

$\displaystyle w(x)=\frac{z(x)}{z(1)}-\left[1+\frac{\beta}{2}(1-x)\right],$
(68)

so that $w(1)=0$. At this stage, we note that it is perfectly equivalent to
look for $w(x)$ or for $z(x)$, and we will see later that many difficulties
can be overcome by using $w(x)$ instead of $z(x)$. The equation for $w(x)$ is
obtained by inserting Eqs. (63) and (68) into Eqs. (65) and (66), in which we
also set $v=1$ and $v^{\prime}=1/y$, with $0<y<1$. Splitting the resulting
integral into a term proportional to $w(x)$ and a second member, we obtain

$\displaystyle\int_{0}^{\infty}\frac{du}{u^{dH}}\left\\{(u+1)^{\beta}w\left(\frac{u+1}{u+y^{-1}}\right)-\frac{u^{\beta}}{G(1)}\left[G\left(\frac{u}{u+1}\right)w\left(\frac{u}{u+y^{-1}}\right)+G\left(\frac{u}{u+y^{-1}}\right)w\left(\frac{u}{u+1}\right)\right]\right\\}=f(y),$
(69)

where the second member $f(y)$ reads

$\displaystyle f(y)=$ $\displaystyle\int_{0}^{\infty}du\ I_{f}(u,y),$ (70)
$\displaystyle I_{f}(u,y)$
$\displaystyle=\frac{(-1)}{u^{dH}}\Big{\\{}(u+1)^{\beta}\left(1+\frac{\beta(1-y)}{2(uy+1)}\right)-\frac{u^{\beta}}{G(1)}\left[G\left(\frac{u}{u+1}\right)\left(1+\frac{\beta}{2(uy+1)}\right)+G\left(\frac{uy}{uy+1}\right)\left(1+\frac{\beta}{2(u+1)}\right)\right]$
$\displaystyle+\frac{u^{\beta}}{G(1)^{2}}G\left(\frac{u}{u+1/y}\right)G\left(\frac{u}{u+1}\right)\left(1+\frac{d}{2}\right)-\frac{d}{2u^{\theta}G(1)}\left[(u+1)^{2H}G\left(\frac{u+1}{u+1/y}\right)-G_{K}(y)\right]\Big{\\}},$
(71)

with the notation

$\displaystyle G_{K}(y)=\begin{cases}G_{s}(y)&\text{if }\theta<1-dH,\\\
G(y)&\text{if }\theta>1-dH.\end{cases}$ (72)

Eq. (69) can be written under the simpler form

$\displaystyle\int_{0}^{1}dx\ K(x,y)\ w(x)=f(y)$ (73)

where the kernel $K(x,y)$ is identified by separating the integrals for the 3
terms in the left-hand side of Eq. (69) and performing a change of variable in
each of them, so that $x$ corresponds to the argument of $w$ (i.e.
$x=(u_{1}+1)/(u_{1}+1/y)$ in the first term, $x=u_{2}/(u_{2}+1/y)$ in the
second one and $x=u_{3}/(u_{3}+1)$ in the third one). This leads to

$\displaystyle
K(x,y)=H(x-y)\frac{du_{1}}{dx}\frac{(u_{1}+1)^{2H-\theta}}{u_{1}^{dH}}-\frac{du_{2}}{dx}\frac{u_{2}^{2H-\theta}}{u_{2}^{dH}G(1)}G\left(\frac{u_{2}}{u_{2}+1}\right)-\frac{du_{3}}{dx}\frac{u_{3}^{2H-\theta}}{u_{3}^{dH}G(1)}G\left(\frac{u_{3}}{u_{3}+1/y}\right),$
(74)

where $H(x)$ is the Heaviside step function, and

$\displaystyle u_{1}=\frac{x-y}{(1-x)y},\hskip
28.45274ptu_{2}=\frac{x}{(1-x)y},\hskip 28.45274ptu_{3}=\frac{x}{1-x}.$ (75)

Technical note for the numerical evaluation of $f(y)$.  Of note, evaluating
numerically $f(y)$ is not straightforward because the integrand $I_{f}(u,y)$
does not always converge rapidly to zero for large $u$ (depending on the value
of $\theta,H,d...$) and a similar problem arises for small $u$. The strategy
to evaluate $f$ numerically consists in identifying analytically the behavior
of $f$ for small/large arguments and using this knowledge in the numerical
program. More precisely, we identify the first terms in the expansion of
$I_{f}(u,y)$ as

$\displaystyle
I_{f}(y,u)\underset{u\to\infty}{\simeq}\frac{a_{1}^{+}}{u^{\nu_{1}^{+}}}+\frac{a_{2}^{+}}{u^{\nu_{2}^{+}}}+\frac{a_{3}^{+}}{u^{\nu_{3}^{+}}}+...,\hskip
28.45274ptI_{f}(y,u)\underset{u\to
0}{\simeq}\frac{a_{1}^{-}}{u^{\nu_{1}^{-}}}+\frac{a_{2}^{-}}{u^{\nu_{2}^{-}}}+\frac{a_{3}^{-}}{u^{\nu_{3}^{-}}}+...$
(76)

Then, $f(y)$ can be evaluated as follows:

$\displaystyle
f(y)=\int_{0}^{1}du\left(I_{f}(y,u)-\sum_{i=1}^{n^{-}}\frac{a_{i}^{-}}{u^{\nu_{i}^{-}}}\right)+\int_{1}^{\infty}du\left(I_{f}(y,u)-\sum_{i=1}^{n^{+}}\frac{a_{i}^{+}}{u^{\nu_{i}^{+}}}\right)+\sum_{i=1}^{n^{-}}\frac{a_{i}^{-}}{1-\nu_{i}^{-}}+\sum_{i=1}^{n^{+}}\frac{a_{i}^{+}}{\nu_{i}^{+}-1},$
(77)

where $n_{+},n^{-}$ are the number of terms for which the expansion of $I_{f}$
is analytically known, and where now the integrals converge rapidly. This
procedure thus require to identify the first order terms of the expansion of
$I_{f}(u,y)$. After some algebra, using the expansion (13) for $G$, we obtain

$\displaystyle\begin{cases}\nu_{1}^{+}=dH+\theta,\hskip
56.9055pt&a_{1}^{+}=d[G_{s}(y)-G_{K}(y)]/[2G(1)]\\\
\nu_{2}^{+}=dH+\theta+2H,\hskip
56.9055pt&a_{2}^{+}=-\kappa_{s}^{2}\left(1+d/2\right)/G(1)^{2}y^{2H}\\\
\nu_{3}^{+}=dH+\theta+1-2H,\hskip 56.9055pt&a_{3}^{+}=0\\\
\nu_{4}^{+}=dH+\theta+1,\hskip
56.9055pt&a_{4}^{+}=\frac{\kappa_{s}}{G(1)}\left(\frac{1}{y^{2H}}+\frac{1}{y}\right)\left[\left(1+\frac{d}{2}\right)H-\frac{\beta}{2}\right]-\frac{d}{G(1)}H\kappa_{s}\left(1/y-1\right)^{2H}\end{cases}$

and, for small arguments:

$\displaystyle\begin{cases}\nu_{1}^{-}=dH+\theta,\hskip
56.9055pt&a_{1}^{-}=d[G(y)-G_{K}(y)]/[2G(1)]\\\ \nu_{2}^{-}=dH,\hskip
56.9055pt&a_{2}^{-}=-\left(1+\beta(1-y)/2\right)\\\
\nu_{3}^{-}=dH+\theta-2H,\hskip
56.9055pt&a_{3}^{-}=2\frac{G(0)}{G(1)}(1+\beta/2)-\frac{G(0)^{2}}{G(1)^{2}}\left(1+\frac{d}{2}\right)\\\
\nu_{4}^{-}=dH-2H+\theta-2\mu,\hskip
56.9055pt&a_{4}^{-}=-\frac{\alpha_{1}^{2}}{G(1)^{2}}\left(1+\frac{d}{2}\right)y^{\mu}\\\
\nu_{5}^{-}=dH-2H+\theta-\mu,\hskip
56.9055pt&a_{5}^{-}=\frac{\alpha_{1}}{G(1)}\left(1+\frac{\beta}{2}\right)(1+y^{\mu})-\frac{G_{0}\alpha_{1}}{G(1)^{2}}(1+y^{\mu})\left(1+\frac{d}{2}\right)\end{cases}$

where $G_{0},\alpha_{1},\mu$ characterize the behavior of $G$ near zero:

$\displaystyle G(x)\underset{x\to 0}{\simeq}G_{0}+\alpha_{1}x^{\mu}+...,\hskip
85.35826pt\mu=\begin{cases}1-2H&(\text{type I})\\\ 1/2-H&(\text{type II})\\\
\end{cases}$

###  D.5 Behavior of $z$ for small arguments and the criterion to select the
correct value of $\theta$

Numerical integration of $z(x)$ tends to show the presence of divergences for
small values of $x$. Here, we give an argument to identify a possible exponent
$\alpha$ such that

$\displaystyle z(x)\underset{x\to 0}{\sim}\frac{1}{x^{\alpha}}.$ (78)

Let us look for non-trivial values of $\alpha$, i.e. which are not contained
in the expansions of $G(x)$. If one looks at Eq. (69) for large $y^{-1}$, we
will find terms $\propto y^{-\alpha}$ that cannot be compensated by any other
term if $\alpha$ is not trivial. The coefficient for this term must then
vanish, leading to

$\displaystyle\int_{0}^{\infty}\frac{dt}{t^{dH}}\left\\{\frac{(t+\tau)^{2H-\theta}}{(t+\tau)^{\alpha}}-\frac{t^{2H-\theta}}{t^{\alpha}}\times\frac{G(t/(t+\tau))}{G(1)}\right\\}y^{-\alpha}+...=0.$
(79)

The values of $\alpha$ that are authorized are those for which the prefactor
of $y^{-\alpha}$ in the above equation vanish, this leads to

$\displaystyle\alpha=\nu+2H-\theta,$ (80)

where $\nu$ is an exponent that does not depend on $\theta$ and satisfies the
equation

$\displaystyle\int_{0}^{\infty}\frac{dt}{t^{dH}}\left\\{\frac{1}{(t+1)^{\nu}}-\frac{1}{t^{\nu}}\times\frac{G(t/(t+1))}{G(1)}\right\\}=0.$
(81)

When $G=G_{s}$ and $d=1$, a solution is $\nu=1-2H$ guerin2016mean  and
$\theta=1-H$ so that $\alpha=H$. This means that most solutions $z(x)$ will
diverge for small arguments. We postulate that the criterion to select the
correct value of $\theta$ is to look for solutions that do not present this
divergence as $1/x^{\alpha}$.

Procedure for the numerical evaluation of $\theta$.  To evaluate $\theta$, we
have proceeded as follows: for a given “test” value of $\theta$, and a given
choice of mesh points $y_{i}$, we evaluated the numerical solution of Eq.
(73). We observed that the obtained values of $z(x)$ behaved as predicted by
(78),(80) and defined a coefficient $C_{\theta}$ so that $z(x)\simeq
C_{\theta}/x^{\nu+2H-\theta}$. Repeating the procedure for several “attempt”
values of $\theta$ enables us to identify the value of $\theta$ for which
$C_{\theta}=0$; we chose this value as the output of our algorithm. We
controlled that the results do not depend on the choice of the mesh.

##  Appendix E Discussion of the validity of our hypotheses

Here, we discuss the validity of our theory. The first issue with our approach
comes from the selection criterium of the correct value of $\theta$ via the
argument that $z(x)$ (or, equivalently, $w(x)$) does not diverge at small
values of $x$. This criterium seems plausible at least when $2H-\theta>0$ in
Eq. (63), since for this value $\rho(t,t^{\prime})\sim t^{2H-\theta}$ diverges
for large $t$. However, this argument is less plausible when $2H-\theta<0$.
This is supported by Fig. 6, where we have compared the theoretical
predictions of $\theta$ to the simulations for all the models presented in the
main text, it is clearly seen that most deviations between theory and
simulations appear when $\theta-2H>0$. We argue that a minimal criterium for
our theory to work is to require that the condition $2H-\theta>0$ is satisfied
at least for the stationary case ($\sigma=\sigma_{s}$), for which
$\theta=1-dH$. This argument leads to the conclusion that we should restrict
ourselves to the case $H>1/(d+2)$.

Figure 6: Comparison of theoretical values of the persistence exponent versus
simulations for all data shown in the main text.

Next, in Section D.2, we have derived an equation that defines $\sigma_{\pi}$.
However, the relation (57) was assumed to hold, and we also assumed that the
right-hand-side of Eq. (61) tends to zero. Here we discuss these hypotheses.
Since the evaluation of $\Sigma$ requires the knowledge of the covariance
after the first-passage, a quantity which is unknown, it seems necessary to
make some ansatz to evaluate it. We have tried two ansatz:

$\displaystyle\Sigma(t_{1},t_{2};u;\tau)\simeq\begin{cases}\sigma_{\pi}^{*}(t_{1},t_{2}|u)&\text{(Ansatz
A: decoupling approx.)}\\\
\text{Cov}[x(u+\tau+t_{1}),x(u+\tau+t_{2})|x(u+\tau)=0;x(\tau)=0]&\text{(Ansatz
B)}\end{cases}$ (82)

In the decoupling approximation (Ansatz A), one assumes that the covariance in
the future of the FPT does not depend on the value of this FPT, this can be
seen as a “mean-field” approximation whose validity can be evaluated self-
consistently. This approximation also reminds us the independent interval
approximation to calculate the duration of intervals between successive zeros
of a smooth Gaussian process, for which the duration of an interval is assumed
to be independent on the duration of the previous interval. In turn, Ansatz B
consists in replacing the condition that $x=0$ was reached for the first time
at $\tau$ by the simpler condition that $x=0$ is observed at $\tau$ (not
necessarily for the first time). This kind of approximation has been proposed
before in the different context of polymer cyclization Sokolov2003 . It is
motivated by the fact that it leads to rather accurate predictions for the
mean square displacement in the future of the first-passage time, see Fig. 5.

In Ansatz $A$, we need to remember that
$\sigma_{\pi}(t,t^{\prime})\simeq\sigma_{0}(t,t^{\prime})$ for large times, at
leading order. It turns out that, in both ansatzs $A$ and $B$, the property
$\Sigma\to\sigma_{s}$ [Eq. (57)] holds, which is an argument in favor of the
validity of our equations. Next, we note that

$\displaystyle\Sigma(t_{1},t_{2};u;\tau)\underset{\tau\to\infty}{\simeq}\begin{cases}\sigma_{\pi}(t_{1},t_{2}|u)&\text{(Ansatz
A),}\\\ \sigma_{s}(t_{1},t_{2}|u)&\text{(Ansatz B).}\\\ \end{cases}$ (83)

This is enough to evaluate the rhs of Eq. (61) as

$\displaystyle R(t_{1},t_{2},T)\simeq\ $ $\displaystyle
T\int_{0}^{T}du\int_{1-u/T}^{\infty}d\tilde{\tau}f(T\tilde{\tau})p(\mathbf{0},T\tilde{\tau}+u|\mathrm{FPT}=T\tilde{\tau})[\Sigma(t_{1},t_{2},u,T\tilde{\tau})-\sigma_{s}(t_{1},t_{2})]$
(84) $\displaystyle=\ $ $\displaystyle
T^{-\theta}\times\begin{cases}\int_{0}^{\infty}du\int_{1}^{\infty}d\tilde{\tau}\frac{f_{0}}{\tilde{\tau}^{1+\theta}}\frac{[\sigma_{\pi}(t_{1},t_{2}|u)-\sigma_{s}(t_{1},t_{2})]}{[2\pi\sigma_{0}(u,u)]^{d/2}}&\text{(Ansatz
A)}\\\
\int_{0}^{\infty}du\int_{1}^{\infty}d\tilde{\tau}\frac{f_{0}}{\tilde{\tau}^{1+\theta}}\frac{[\sigma_{s}(t_{1},t_{2}|u)-\sigma_{s}(t_{1},t_{2})]}{[2\pi\sigma_{s}(u,u)]^{d/2}}&\text{(Ansatz
B)}\end{cases}$ (85)

In Eq. (84) we have set $\tau=T\tilde{\tau}$ and the next evaluation follows
from taking the limit $T\to\infty$ of all terms. The integrals in the above
equation converge as soon as $\theta>0$, $dH<1$ (condition obtained for small
$u$), and $(d+2)H>1$ (condition obtained for large $u$, for which
$\sigma_{s}(t_{1},t_{2}|u)-\sigma_{s}(t_{1},t_{2})\propto 1/u^{2H}$ from the
property (16), and the same holds for $\sigma_{\pi}$ when one assumes that it
has also long-time stationary increments). This means that, if $(d+2)H>1$, the
leading order of the lhs of Eq. (61) is $1/T^{\theta}$ for large $T$, which
tends to zero, so that our equation for the covariance is correct (at least
from the point of view of these investigations with our two ansatzs). Hence,
we expect that our theory holds when $H>1/(d+2)$, this validity criterium is
the same as the one found above.

In order to investigate the behavior of $R$ for $H<1/(d+2)$, it is necessary
to show that

$\displaystyle\Sigma(t_{1},t_{2},u,\tau)\underset{u,\tau\to\infty}{\simeq}\begin{cases}\sigma_{s}(t_{1},t_{2})-\frac{\kappa_{s}^{2}(t_{1}t_{2})^{2H}}{G(1)u^{2H}}\hskip
28.45274pt&\text{(Ansatz A}),\\\ \sigma_{s}(t_{1},t_{2})-\frac{\kappa_{s}^{2}\
(t_{1}t_{2})^{2H}}{\sigma(u,u|\tau)}\hskip 28.45274pt&\text{(Ansatz
B}),\end{cases}$ (86)

where the result for Ansatz $B$ comes from general formulas on conditional
covariances Eaton1983  and the property (16). If $(d+2)H<1$, we thus evaluate
the rhs of Eq. (61) as

$\displaystyle R(t_{1},t_{2},T)$
$\displaystyle=T^{2}\int_{1-\tilde{u}}^{\infty}d\tilde{\tau}f(T\tilde{\tau})p(\mathbf{0},T\tilde{\tau}+T\tilde{u}|\mathrm{FPT}=T\tilde{\tau})[\Sigma(t_{1},t_{2},T\tilde{u},T\tilde{\tau})-\sigma_{s}(t_{1},t_{2})]$
(87)
$\displaystyle\simeq-\frac{T^{1-(d+2)H-\theta}}{(2\pi)^{d/2}}\kappa_{s}^{2}\times\begin{cases}\int_{0}^{1}d\tilde{u}\int_{1-\tilde{u}}^{\infty}d\tilde{\tau}\frac{f_{0}}{\tilde{\tau}^{1+\theta}}\frac{1}{[G(1)u]^{(d+2)H}}&\text{(Ansatz
A)},\\\
\int_{0}^{1}d\tilde{u}\int_{1-\tilde{u}}^{\infty}d\tilde{\tau}\frac{f_{0}}{\tilde{\tau}^{1+\theta}}\frac{1}{[\sigma(u,u|\tau)]^{{d+2}H}}&(\text{Ansatz
B}).\end{cases}$ (88)

These integrals exist when $0<\theta<1$ and $(d+2)H<1$ (which is assumed
here), so that in this case $R\propto T^{1-(d+2)H-\theta}$. This means that
when $\theta<1-(d+2)H$, our equation for the covariance is wrong since the rhs
of (61) does not vanish for large $T$. Since there is no reason why $\theta$
should be larger than $1-(d+2)H$, we consider that using our formalism for
$(d+2)H<1$ is not advised.

In Fig. 7, we present simulation results for a value of $H$ which does not
satisfy this criterium ($H=1/4,d=1$). It is seen that the curve $\theta(T)$ is
qualitatively correct, but there is a clear difference between the predicted
value of $\theta$ and the measured one. However, from the arguments above, we
are unsure of our criterium as soon as $\theta>2H=1/2$, and we know that our
equation for the covariance is wrong for $\theta<1-3H=0.25$. This leaves a
very small range of values of $\theta$ in which the theory could be expected
to be correct (which does not even include the stationary state $T=1$) and
this is probably the reason of the differences between theory and simulations.

Figure 7: Results for $H=1/4$ and $d=1$ (which lies outside the validity
criterium $H>1/(d+2)$ of our theory, see text).

##  Appendix F Perturbative analysis of the formalism [Derivation of Eqs.
(7),(8)]

Here, we perform an analysis at perturbative order of our equations in
dimension $d=1$. We assume $H=1/2+\varepsilon$ with $\varepsilon$ small and we
start directly with the equation (66) for $\rho(u,v)$. We assume $\theta>1-H$
(but we expect that the expansion in powers of $\varepsilon$ will be the same
for $\theta>1-H$ or $\theta<1-H$):

$\displaystyle\int_{0}^{\infty}\frac{du}{u^{H}}\Bigg{\\{}\rho(u+v,u+v^{\prime})-$
$\displaystyle\frac{\sigma(u+v,u)\rho(u+v^{\prime},u)}{\sigma(u,u)}-\frac{\sigma(u+v^{\prime},u)\rho(u+v,u)}{\sigma(u,u)}+\rho(u,u)\frac{\sigma(u+v,u)\sigma(u+v^{\prime},u)}{\sigma(u,u)^{2}}\Bigg{\\}}$
$\displaystyle=\int_{0}^{\infty}\frac{du}{2u^{H+\theta}}\left[\sigma(u+v,u+v^{\prime})-\frac{\sigma(u+v,u)\sigma(u+v^{\prime},u)}{\sigma(u,u)}-\sigma(v,v^{\prime})\right]\equiv
C(v,v^{\prime}),$ (89)

where $C$ is defined by this equation. Following our formalism, we will look
for symmetric solutions with $\rho(u,u)=u^{2H-\theta}$ and that our criterium
to select the correct value of $\theta$ is that $\rho(u,v)$ does not diverge
with $v$ as $v\to 0$. Let us write the expansion of $\theta$ as:

$\displaystyle\theta=1-H+\varepsilon\delta_{1}+\varepsilon^{2}\delta_{2}+...,$
(90)

with $\delta_{1}\varepsilon>0$ to have $\theta>1-H$. Now, let us introduce the
following notations

$\displaystyle\sigma_{s}(t,t^{\prime})=\text{min}(t,t^{\prime})+\varepsilon\sigma_{s,1}(t,t^{\prime})+\varepsilon^{2}\sigma_{s,2}(t,t^{\prime})+...$
$\displaystyle\sigma(t,t^{\prime})=\sigma_{s}(t,t^{\prime})+\varepsilon\omega_{1}(t,t^{\prime})+\varepsilon^{2}\omega_{2}(t,t^{\prime})+...$
(91)
$\displaystyle\rho(u,v)=\rho_{0}(u,v)+\varepsilon\rho_{1}(u,v)+\varepsilon^{2}\rho_{2}(u,v)...$
(92)

Here, $\omega_{i}$ represents the deviation of $\sigma$ with respect to
$\sigma_{s}$ at order $i$. Note that the leading order term for $\rho$ is of
order zero because of the normalization
$\rho(u,u)=u^{2H-\theta}\simeq\sqrt{u}$ when $H\to 1/2$.

###  F.1 First order

To estimate the limit of $C$ for $\varepsilon\to 0$, let us introduce a
parameter $L$ (finite but large) and consider the calculation

$\displaystyle C(v,v^{\prime})$
$\displaystyle=\int_{0}^{L}\frac{du}{2u^{H+\theta}}\left[\sigma(u+v,u+v^{\prime})-\frac{\sigma(u+v,u)\sigma(u+v^{\prime},u)}{\sigma(u,u)}-\sigma(v,v^{\prime})\right]+\int_{L}^{\infty}\frac{du}{2u^{H+\theta}}\left[\sigma_{s}(v,v^{\prime})-\sigma(v,v^{\prime})\right],$
$\displaystyle=\int_{0}^{L}\frac{du}{2u^{H+\theta}}\left[\sigma(u+v,u+v^{\prime})-\frac{\sigma(u+v,u)\sigma(u+v^{\prime},u)}{\sigma(u,u)}-\sigma(v,v^{\prime})\right]+\frac{L^{-\varepsilon\delta_{1}}}{2\varepsilon\delta_{1}}\left[\sigma_{s}(v,v^{\prime})-\sigma(v,v^{\prime})\right],$
$\displaystyle\underset{\varepsilon\to
0}{\simeq}-\frac{\omega_{1}(v,v^{\prime})}{2\delta_{1}}.$ (93)

We thus obtain

$\displaystyle\int_{0}^{\infty}\frac{du}{u^{1/2}}\Bigg{\\{}\rho_{0}(u+v,u+v^{\prime})-\rho_{0}(u+v^{\prime},u)-\rho_{0}(u,u+v^{\prime})+\rho_{0}(u,u)\Bigg{\\}}=-\frac{\omega_{1}(v,v^{\prime})}{2\delta_{1}}.$
(94)

Let us take the derivative with respect to $v$ and $v^{\prime}$ of this
equation:

$\displaystyle\int_{0}^{\infty}\frac{du}{u^{1/2}}\rho_{0}^{(1,1)}(u+v,u+v^{\prime})=-\frac{\omega_{1}^{(1,1)}(v,v^{\prime})}{2\delta_{1}},$
(95)

with the notation $G^{(1,1)}=\partial_{v}\partial_{v^{\prime}}G$ for all
functions $G(v,v^{\prime})$. The solution of this equation can be found by the
method of superposition of solutions. Consider the simpler equation with a
single complex exponential in the second member

$\displaystyle\int_{0}^{\infty}\frac{du}{u^{1/2}}\rho_{0}^{(1,1)}(u+v,u+v^{\prime})=e^{i(\omega
v+\omega^{\prime}v^{\prime})},$ (96)

for which an obvious solution is

$\displaystyle\rho_{0}^{(1,1)}(v,v^{\prime})=e^{i(\omega
v+\omega^{\prime}v^{\prime})}\frac{1}{\int_{0}^{\infty}\frac{du}{u^{1/2}}e^{i(\omega+\omega^{\prime})u}}.$
(97)

Introducing the double Fourier transform and its inverse for all functions
$f(t,t^{\prime})$ as

$\displaystyle\hat{f}(\omega,\omega^{\prime})=\int_{-\infty}^{\infty}dt\int_{-\infty}^{\infty}dt^{\prime}e^{-i(\omega
t+\omega^{\prime}t^{\prime})}f(t,t^{\prime}),\hskip
28.45274ptf(t,t^{\prime})=\frac{1}{4\pi^{2}}\int_{-\infty}^{\infty}d\omega\int_{-\infty}^{\infty}d\omega^{\prime}e^{i(\omega
t+\omega^{\prime}t^{\prime})}\hat{f}(\omega,\omega^{\prime}),$ (98)

we see that

$\displaystyle\rho_{0}^{(1,1)}(v,v^{\prime})=-\frac{1}{8\pi^{2}\delta_{1}}\int_{-\infty}^{\infty}d\omega\int_{-\infty}^{\infty}d\omega^{\prime}\
\widehat{\omega_{1}^{(1,1)}}(\omega,\omega^{\prime})e^{i(\omega
v+\omega^{\prime}v^{\prime})}\frac{1}{\int_{0}^{\infty}\frac{du}{u^{1/2}}e^{i(\omega+\omega^{\prime})u}}$
(99)

is a solution of (95) constructed by superposition of solutions. Integrating
with respect to $v,v^{\prime}$ and using
$\rho_{1}(v,0)=\rho_{1}(0,v^{\prime})=0$, we obtain

$\displaystyle\rho_{0}(\xi,\xi^{\prime})=\int_{0}^{\xi}dv\int_{0}^{\xi^{\prime}}dv^{\prime}\rho_{0}^{(1,1)}(v,v^{\prime})$
$\displaystyle=-\frac{1}{8\pi^{2}\delta_{1}}\int_{-\infty}^{\infty}d\omega\int_{-\infty}^{\infty}d\omega^{\prime}\int_{0}^{\infty}dt\int_{0}^{\infty}dt^{\prime}\omega_{1}^{(1,1)}(t,t^{\prime})\frac{(e^{i\omega(\xi-t)}-e^{-i\omega
t})(e^{i\omega^{\prime}(\xi^{\prime}-t^{\prime})}-e^{-i\omega^{\prime}t^{\prime}})}{(-\omega\omega^{\prime})}\frac{\sqrt{-i(\omega+\omega^{\prime})}}{\sqrt{\pi}}.$
(100)

We introduce a function $K_{0}(t,t^{\prime})$ and its double Fourier transform

$\displaystyle
K_{0}(t,t^{\prime})=\left(\frac{H(t^{\prime}-t)}{\pi\sqrt{t}}+\frac{H(t-t^{\prime})}{\pi\sqrt{t^{\prime}}}\right)H(t)H(t^{\prime}),\hskip
28.45274pt\widehat{K}_{0}(\omega,\omega^{\prime})=-\frac{\sqrt{i(\omega+\omega^{\prime})}}{\omega\omega^{\prime}\sqrt{\pi}},$
(101)

with $H$ the Heaviside step function. Therefore, if we change
$\omega\to-\omega$ and $\omega^{\prime}\to-\omega^{\prime}$ in (100) and
integrate over $\omega,\omega^{\prime}$ before $t,t^{\prime}$, we recognize
the double inverse Fourier transform of $K$ and we obtain

$\displaystyle\rho_{0}(\xi,\xi^{\prime})$
$\displaystyle=-\frac{1}{2\delta_{1}}\int_{0}^{\infty}dt\int_{0}^{\infty}dt^{\prime}\omega_{1}^{(1,1)}(t,t^{\prime})[K_{0}(t-\xi,t^{\prime}-\xi^{\prime})-K_{0}(t-\xi,t^{\prime})-K_{0}(t,t^{\prime}-\xi^{\prime})+K_{0}(t,t^{\prime})].$
(102)

Now, using (101) this equation can be integrated once guerin2016mean , leading
to

$\displaystyle\rho_{0}(\xi,\xi^{\prime})=\frac{1}{2\pi\delta_{1}}\int_{0}^{\infty}\frac{dt}{\sqrt{t}}[\omega_{1}^{*}(t+\xi,t+\xi^{\prime})-\omega_{1}^{*}(t+\xi,t)-\omega_{1}^{*}(t,t+\xi^{\prime})+\omega_{1}^{*}(t,t)],$
(103)

with $\omega_{1}^{*}$ is a symetrized first derivative of $\omega_{1}$, i.e.

$\displaystyle\omega_{1}^{*}(t,t^{\prime})=(\partial_{t}+\partial_{t^{\prime}})\omega_{1}(t,t^{\prime}).$
(104)

Eq. (103) is the solution of Eq. (94). For the quenched fBm, we can calculate

$\displaystyle\omega_{1}^{*}(t,t^{\prime})=\lambda[2\ln(t+t^{\prime})-\ln
t-\ln t^{\prime}],$ (105)

with $\lambda=1-T$. Inserting this value into Eq. (103) and performing the
integral, we obtain

$\displaystyle\rho_{0}(\xi,\xi^{\prime})=\frac{\lambda\sqrt{2}}{\delta_{1}}[\sqrt{\xi+\xi^{\prime}}-\sqrt{\xi}-\sqrt{\xi^{\prime}}].$
(106)

Comparing the above equation with the property $\rho_{0}(\xi,\xi)=\sqrt{\xi}$
(which holds by the very definiton of $\rho_{0}$) enables us to identify
$\delta_{1}$:

$\displaystyle\delta_{1}=-2\lambda(\sqrt{2}-1).$ (107)

This is consistent with (and generalizes) the results of Ref Krug1997  for
$\lambda=1$ and $\lambda=0$. It is also instructive to give the value of
$z(x)$ at this order:

$\displaystyle
z(x)=-\frac{\sqrt{1+x}-\sqrt{1}-\sqrt{x}}{(2-\sqrt{2})\sqrt{x}}.$ (108)

###  F.2 Second order (quenched fBm)

At second order, we need the next-to-leading order of the expansion of the
term $C$ in Eq. (89):

$\displaystyle
C(v,v^{\prime})=\int_{0}^{L}\frac{du}{2u^{1+\varepsilon\delta_{1}}}\left[\sigma(u+v,u+v^{\prime})-\frac{\sigma(u+v,u)\sigma(u+v^{\prime},u)}{\sigma(u,u)}-\sigma(v,v^{\prime})\right]+\frac{L^{-\varepsilon\delta_{1}-\varepsilon^{2}\delta_{2}}}{2[\varepsilon\delta_{1}+\varepsilon^{2}\delta_{2}]}\left[\sigma_{s}(v,v^{\prime})-\sigma(v,v^{\prime})\right]$
$\displaystyle\simeq\varepsilon\int_{0}^{L}\frac{du}{2u}\left[\sigma_{1}(u+v,u+v^{\prime})-\sigma_{1}(u+v,u)-\sigma_{1}(u,u+v^{\prime})+\sigma_{1}(u,u)-\sigma_{1}(v,v^{\prime})\right]+\frac{1-\varepsilon\delta_{1}\ln
L}{2\varepsilon\delta_{1}(1+\varepsilon\delta_{2}/\delta_{1})}\left[-\varepsilon\omega_{1}-\varepsilon^{2}\omega_{2}\right]$
$\displaystyle\simeq-\frac{\omega_{1}(v,v^{\prime})}{2\delta_{1}}+\varepsilon
C_{2}$ (109)

with

$\displaystyle C_{2}(v,v^{\prime})=\int_{0}^{L}\frac{du}{2u}$
$\displaystyle\left[\omega_{1}(u+v,u+v^{\prime})-\omega_{1}(u+v,u)-\omega_{1}(u,u+v^{\prime})+\omega_{1}(u,u)-\omega_{1}(v,v^{\prime})\right]$
$\displaystyle+\frac{\omega_{1}(v,v^{\prime})\ln
L}{2}+\omega_{1}(v,v^{\prime})\frac{\delta_{2}}{2\delta_{1}^{2}}-\frac{\omega_{2}(v,v^{\prime})}{2\delta_{1}}$
(110)

Using $\ln L=\int_{1}^{L}du/u$ in the above expression, we see that $C_{2}$ is
independent of $L$ and reads

$\displaystyle C_{2}(v,v^{\prime})=$
$\displaystyle\int_{0}^{1}\frac{du}{2u}\left[\omega_{1}(u+v,u+v^{\prime})-\omega_{1}(u+v,u)-\omega_{1}(u,u+v^{\prime})+\omega_{1}(u,u)-\omega_{1}(v,v^{\prime})\right]+$
$\displaystyle\int_{1}^{\infty}\frac{du}{2u}\left[\omega_{1}(u+v,u+v^{\prime})-\omega_{1}(u+v,u)-\omega_{1}(u,u+v^{\prime})+\omega_{1}(u,u)\right]+\omega_{1}(v,v^{\prime})\frac{\delta_{2}}{2\delta_{1}^{2}}-\frac{\omega_{2}(v,v^{\prime})}{2\delta_{1}}.$
(111)

Expanding Eq. (89) at second order, we find that the equation for $\rho_{1}$
is

$\displaystyle\int_{0}^{\infty}\frac{du}{u^{1/2}}\Bigg{\\{}\rho_{1}(u+v,u+v^{\prime})-\rho_{1}(u+v,u)-\rho_{1}(u,u+v^{\prime})+\rho_{1}(u,u)\Bigg{\\}}=A_{1}+A_{2}+C_{2}\equiv
A(v,v^{\prime})$ (112) $\displaystyle
A_{1}=\int_{0}^{\infty}\frac{du}{u^{3/2}}\Big{\\{}[\sigma_{1}(u+v,u)-\sigma_{1}(u,u)][\rho_{0}(u+v^{\prime},u)-\rho_{0}(u,u)]+[\sigma_{1}(u+v^{\prime},u)-\sigma_{1}(u,u)][\rho_{0}(u+v,u)-\rho_{0}(u,u)]\Big{\\}}$
$\displaystyle A_{2}=\int_{0}^{\infty}du\left(\frac{\ln
u}{u^{1/2}}\right)\Bigg{\\{}\rho_{0}(u+v,u+v^{\prime})-\rho_{0}(u+v,u)-\rho_{0}(u,u+v^{\prime})+\rho_{0}(u,u)\Bigg{\\}}$
(113)

Eq. (112) is the same integral equation that appeared at first order, its
solution with no divergence at the origin is

$\displaystyle\rho_{1}(\xi,\xi)=-\frac{1}{\pi}\int_{0}^{\infty}\frac{dt}{\sqrt{t}}[A^{*}(t+\xi,t+\xi)-2A^{*}(t+\xi,t)+A^{*}(t,t)].$
(114)

Performing the integrations, we obtain

$\displaystyle\rho_{1}(\xi,\xi)=\sqrt{\xi}[\alpha+\Omega\ln\xi],\hskip
28.45274pt\Omega=3+2(-1+\sqrt{2})\lambda=3-\delta_{1}.$ (115)
$\displaystyle\alpha=-4+c+\frac{(1+\sqrt{2})\delta_{2}}{2\lambda}+\lambda(4-4\sqrt{2}+b-6\ln
2+\sqrt{2}\ln 4)-\sqrt{2}\ln 8$ (116)

where the numerical constants $b$ and $c$ are defined by the integrals

$\displaystyle c=$
$\displaystyle-\frac{(2+\sqrt{2})}{\pi}\int_{0}^{\infty}\frac{d\tau}{\sqrt{\tau}}\int_{0}^{\infty}\frac{du}{u^{3/2}}\frac{\partial}{\partial\tau}\Big{\\{}\left(\sqrt{\tau+u}-\sqrt{1+\tau+u}-\sqrt{\tau+2u}+\sqrt{1+\tau+2u}\
\right)\times$
$\displaystyle\left[(u+\tau)\ln(u+\tau)-\tau\ln\tau+(1+\tau)\ln(1+\tau)-(1+u+\tau)\ln(1+u+\tau)\right]\Big{\\}}$
(117) $\displaystyle b=$
$\displaystyle-\frac{2}{\pi(-2+\sqrt{2})}\int_{0}^{\infty}\frac{d\tau}{\sqrt{\tau}}\int_{0}^{\infty}\frac{du}{u^{3/2}}\frac{\partial}{\partial\tau}\Big{\\{}\left(\sqrt{\tau+u}-\sqrt{1+\tau+u}-\sqrt{\tau+2u}+\sqrt{1+\tau+2u}\
\right)\times$
$\displaystyle\left[(u+\tau)\ln(u+\tau)-(1+u+\tau)\ln(1+u+\tau)-(2u+\tau)\ln(2u+\tau)+(1+2u+\tau)\ln(1+2u+\tau)\right]\Big{\\}}$
(118)

To comment on the value of $\Omega$, let us remind that

$\displaystyle\rho(\xi,\xi)=\xi^{\beta}=\xi^{2H-\theta}=\xi^{3H-1-\varepsilon\delta_{1}+...}=\sqrt{\xi}+\varepsilon(3-\delta_{1})\sqrt{\xi}\ln\xi+\mathcal{O}(\varepsilon^{2}).$
(119)

If our theory is consistent the coefficient of $\varepsilon\sqrt{\xi}\ln\xi$
in the above equation is set by $\delta_{1}$, whatever the value of
$\delta_{2}$. The fact that we find the same value with the above argument
(119) or by the complete calculation (115) argues in favor of the consistency
of our approach. We also see that there is no term $\varepsilon\sqrt{\xi}$ at
first order in the above expression (119), which imposes to set $\alpha=0$;
this equality can be satisfied by adjusting the value of $\delta_{2}$. The
result is

$\displaystyle\delta_{2}=a_{1}\lambda(a_{2}+\lambda),$ (120)

with the constants

$\displaystyle a_{1}=-2(-1+\sqrt{2})(4-4\sqrt{2}+b-6\ln 2+\sqrt{2}\ln 4)\simeq
1.77,\hskip 28.45274pta_{2}=\frac{-4+c-\sqrt{2}\ln 8}{4-4\sqrt{2}+b-6\ln
2+\sqrt{2}\ln 4}\simeq 1.28,$ (121)

where the numerical values necessitate the numerical evaluation of the
integrals in Eqs. (117) and (118). To summarize, our result at second order
for the quenched fBm is

$\displaystyle\theta=1-H+\varepsilon\delta_{1}+\varepsilon^{2}\delta_{2},\hskip
28.45274pt \delta_{1}=-2\lambda(\sqrt{2}-1)\hskip
28.45274pt\delta_{2}=1.77\lambda(1.28+\lambda)\hskip 28.45274pt\lambda=1-T.$
(122)

###  F.3 Second order for the fBm constrained on its past

For the fBm constrained on its past, we have again to evaluate the second
member $C$ in Eq. (89), which reads:

$\displaystyle C(v,v^{\prime})$
$\displaystyle=\int_{0}^{L}\frac{du}{2u^{H+\theta}}\left[\sigma(u+v,u+v^{\prime})-\frac{\sigma(u+v,u)\sigma(u+v^{\prime},u)}{\sigma(u,u)}-\sigma(v,v^{\prime})\right]+\int_{L}^{\infty}\frac{du}{2u^{H+\theta}}\left[\sigma_{s}(v,v^{\prime})-\sigma(v,v^{\prime})\right].$
(123)

Now, since $\sigma=\sigma_{s}$ at first order, and
$\theta=1-H+\varepsilon^{2}\delta_{2}$, we get:

$\displaystyle C(v,v^{\prime})$
$\displaystyle=\int_{0}^{L}\frac{du}{2u^{H+\theta}}\left[\sigma(u+v,u+v^{\prime})-\frac{\sigma(u+v,u)\sigma(u+v^{\prime},u)}{\sigma(u,u)}-\sigma(v,v^{\prime})\right]+\frac{L^{-\varepsilon^{2}\delta_{2}}}{2\varepsilon^{2}\delta_{2}}\left[\sigma_{s}(v,v^{\prime})-\sigma(v,v^{\prime})\right]$
$\displaystyle\simeq-\frac{\omega_{2}(v,v^{\prime})}{2\delta_{2}}.$ (124)

This means that the equation for $\rho_{0}$ is essentially the same as in the
first order case, Eq. (94) as soon as one replaces $\omega_{1}/\delta_{1}$ by
$\omega_{2}/\delta_{2}$. From Eq. (103) we deduce the solution for $\rho_{0}$

$\displaystyle\rho_{0}(\xi,\xi^{\prime})=\frac{1}{2\pi\delta_{2}}\int_{0}^{\infty}\frac{dt}{\sqrt{t}}[\omega_{2}^{*}(t+\xi,t+\xi^{\prime})-\omega_{2}^{*}(t+\xi,t)-\omega_{2}^{*}(t,t+\xi^{\prime})+\omega_{2}^{*}(t,t)].$
(125)

We can identify $\omega_{2}$ by expanding (27) at second order

$\displaystyle\omega_{2}(t,t^{\prime})=-\int_{0}^{\infty}dx\frac{(t+xt^{\prime})\ln(t+xt^{\prime})-t\ln
t-(t^{\prime}x)\ln(t^{\prime}x)}{x(1+x)}.$ (126)

The integral can be performed after taking the derivatives with respect to $t$
and $t^{\prime}$

$\displaystyle\omega_{2}^{*}(t,t^{\prime})=(\partial_{t}+\partial_{t^{\prime}})\omega_{2}=-\frac{1}{6}(2\pi^{2}+3\ln^{2}(t^{\prime}/t)).$
(127)

Inserting this value into Eq. (125) leads to

$\displaystyle\rho_{0}(\xi,\xi^{\prime})=\frac{2}{\delta_{2}}\left[(\sqrt{\xi}+\sqrt{\xi^{\prime}})\ln(\sqrt{\xi}+\sqrt{\xi^{\prime}})-\sqrt{\xi}\ln\sqrt{\xi}-\sqrt{\xi^{\prime}}\ln\sqrt{\xi^{\prime}}\right]=\frac{1}{4\pi\delta_{2}}\sqrt{\xi}\
z(\xi/\xi^{\prime})$ (128)

with

$\displaystyle
z(x)=\frac{\pi}{\sqrt{x}}[8(1+\sqrt{x})\ln(1+\sqrt{x})-4\sqrt{x}\ln x]$ (129)

Strikingly we observe a logarithmic divergence of $z(x)$ around $x=0$, as
$z(x)\propto 1-2\ln x$. Again, the argument $\rho_{0}(\xi,\xi)=\sqrt{\xi}$
enables us to identify the value of $\delta_{2}=4\ln 2$, and thus the result
of our perturbation theory for the fBm constrained on its past reads

$\displaystyle\theta\simeq 1-H+(4\ln 2)\varepsilon^{2}+...$ (130)

##  Appendix G Result in the limit of high initial temperature, $T=\infty$

Let us consider the quenched fBm in the limit $T\to\infty$, which corresponds
physically to the case that the initial temperature is much larger than the
temperature of the dynamics. Thus, this corresponds to models starting with a
fully “disordered” state. We see that the limiting value of the scaling
function becomes

$\displaystyle G_{\infty}(u)\equiv
G_{T=\infty}(u)=\frac{-(1+u)^{2H}+(1+u^{2H})}{[-2^{2H}+2]u^{2H}}.$ (131)

thus we see that the behavior of $G_{T=\infty}(u)$ near $u=1$ is regular.
Hence, according to Eq. (18), after the Lamperti transform we obtain a smooth
Gaussian process whose density of zero crossing can be analyzed with the
independent interval approximation. Denoting again by $a(T)$ the correlator of
the process after Lamperti transform,
$a_{\infty}(T)=e^{-HT}G_{\infty}(e^{-HT})$, the result of the IIA is that
$\theta$ is the smallest positive number that satisfies the equation
$F(-\theta)=0$, where

$\displaystyle F(s)=1+\frac{\pi\
s}{2|a_{\infty}^{\prime\prime}(0)|}\left[1-\frac{2s}{\pi}\int_{0}^{\infty}dTe^{-sT}\text{asin}(a_{\infty}(T))\right].$
(132)

The value of $\theta$ deduced from this approximation for $T=\infty$ are close
to $\theta\simeq 0.11\pm 0.01$ for all values of $H$. This means that the
perturbation theory around $H=1/2$ will fail at high temperatures, since the
weakly non-Markovian limit and the limit of large initial temperature cannot
be inverted.

## References

  * (1) Majumdar, S. N.  Persistence in nonequilibrium systems.  _Current Science_ 77 (1999).
  * (2) Bray, A. J., Majumdar, S. N. & Schehr, G.  Persistence and first-passage properties in nonequilibrium systems.  _Adv. Phys._ 62, 225–361 (2013).
  * (3) Wong, G. P., Mair, R. W., Walsworth, R. L. & Cory, D. G.  Measurement of persistence in 1d diffusion.  _Phys. Rev. Lett._ 86, 4156 (2001).
  * (4) Tam, W., Zeitak, R., Szeto, K. & Stavans, J.  First-passage exponent in two-dimensional soap froth.  _Phys. Rev. Lett._ 78, 1588 (1997).
  * (5) Dougherty, D. _et al._ Experimental persistence probability for fluctuating steps.  _Phys. Rev. Lett._ 89, 136102 (2002).
  * (6) Merikoski, J., Maunuksela, J., Myllys, M., Timonen, J. & Alava, M. J.  Temporal and spatial persistence of combustion fronts in paper.  _Phys. Rev. Lett._ 90, 024501 (2003).
  * (7) Soriano, J., Braslavsky, I., Xu, D., Krichevsky, O. & Stavans, J.  Universality of persistence exponents in two-dimensional ostwald ripening.  _Phys. Rev. Lett._ 103, 226101 (2009).
  * (8) Majumdar, S., Sire, C., Bray, A. & Cornell, S. J.  Nontrivial exponent for simple diffusion.  _Phys. Rev. Let.t_ 77, 2867–2870 (1996).
  * (9) Derrida, Hakim, V. & Zeitak.  Persistent spins in the linear diffusion approximation of phase ordering and zeros of stationary gaussian processes.  _Phys. Rev. Lett._ 77, 2871–2874 (1996).
  * (10) Derrida, B., Hakim, V. & Pasquier, V.  Exact first-passage exponents of 1d domain growth: relation to a reaction-diffusion model.  _Phys. Rev. Lett._ 75, 751 (1995).
  * (11) Majumdar, S. & Sire, C.  Survival probability of a gaussian non-markovian process: Application to the t=0 dynamics of the ising model.  _Phys. Rev. Lett._ 77, 1420–1423 (1996).
  * (12) Majumdar, S. N. & Bray, A. J.  Spatial persistence of fluctuating interfaces.  _Phys. Rev. Lett._ 86, 3700–3 (2001).
  * (13) Majumdar, S., Bray, A., Cornell, S. & Sire, C.  Global persistence exponent for nonequilibrium critical dynamics.  _Phys. Rev. Lett._ 77, 3704 (1996).
  * (14) Ehrhardt, G. C. & Bray, A. J.  Series expansion calculation of persistence exponents.  _Phys. Rev. Lett._ 88, 070601 (2002).
  * (15) Lee, B. P. & Rutenberg, A. D.  Persistence, poisoning, and autocorrelations in dilute coarsening.  _Phys. Rev. Lett._ 79, 4842 (1997).
  * (16) Constantin, M. _et al._ Infinite family of persistence exponents for interface fluctuations.  _Phys. Rev. Lett._ 91, 086103 (2003).
  * (17) Krug, J. _et al._ Persistence exponents for fluctuating interfaces.  _Phys. Rev. E_ 56, 2702–2712 (1997).
  * (18) Molchan, G.  Maximum of a fractional brownian motion: Probabilities of small values.  _Commun. Math. Phys._ 205, 97–111 (1999).
  * (19) Derrida, B., Bray, A. & Godreche, C.  Non-trivial exponents in the zero temperature dynamics of the 1d ising and potts models.  _J. Phys. A: Math. Gen._ 27, L357 (1994).
  * (20) Redner, S.  _A guide to First- Passage Processes_ (Cambridge University Press, Cambridge, England, 2001).
  * (21) Condamin, S., Bénichou, O., Tejedor, V., Voituriez, R. & Klafter, J.  First-passage times in complex scale-invariant media.  _Nature_ 450, 77–80 (2007).
  * (22) Bénichou, O., Chevalier, C., Klafter, J., Meyer, B. & Voituriez, R.  Geometry-controlled kinetics.  _Nat. Chem._ 2, 472–477 (2010).
  * (23) Schuss, Z., Singer, A. & Holcman, D.  The narrow escape problem for diffusion in cellular microdomains.  _Proc. Natl. Acad. Sci. U. S. A._ 104, 16098–103 (2007).
  * (24) Godec, A. & Metzler, R.  Universal proximity effect in target search kinetics in the few-encounter limit.  _Phys. Rev. X_ 6, 041037 (2016).
  * (25) Delorme, M. & Wiese, K. J.  Maximum of a fractional brownian motion: analytic results from perturbation theory.  _Phys. Rev. Lett._ 115, 210601 (2015).
  * (26) Guérin, T., Levernier, N., Bénichou, O. & Voituriez, R.  Mean first-passage times of non-markovian random walkers in confinement.  _Nature_ 534, 356–359 (2016).
  * (27) Levernier, N., Dolgushev, M., Bénichou, O., Voituriez, R. & Guérin, T.  Survival probability of stochastic processes beyond persistence exponents.  _Nat. Comm._ 10, 1–7 (2019).
  * (28) Dolgushev, M., Guérin, T., Blumen, A., Bénichou, O. & Voituriez, R.  Contact kinetics in fractal macromolecules.  _Phys. Rev. Lett._ 115, 208301 (2015).
  * (29) Guérin, T., Bénichou, O. & Voituriez, R.  Non-markovian polymer reaction kinetics.  _Nat. chem._ 4, 568–573 (2012).
  * (30) Levernier, N., Bénichou, O., Guérin, T. & Voituriez, R.  Universal first-passage statistics in aging media.  _Phys. Rev. E_ 98, 022125 (2018).
  * (31) Burkhardt, T.  Semiflexible polymer in the half-plane and statistics of the integral of a brownian curve.  _J. Phys. A: Math. Gen._ 26, L1157–L1162 (1993).
  * (32) De Smedt, G., Godreche, C. & Luck, J.  Partial survival and inelastic collapse for a randomly accelerated particle.  _Europhys. Lett._ 53, 438 (2001).
  * (33) Poplavskyi, M. & Schehr, G.  Exact persistence exponent for the 2 d-diffusion equation and related kac polynomials.  _Phys. Rev. Lett._ 121, 150601 (2018).
  * (34) Dornic, I.  Universal painlevé vi probability distribution in pfaffian persistence and gaussian first-passage problems with a sech-kernel.  _arXiv preprint arXiv:1810.06957_ (2018).
  * (35) Watson, A.  Persistence pays off in defining history of diffusion.  _Science_ 274, 919–920 (1996).
  * (36) Newman, T. & Loinaz, W.  Critical dimensions of the diffusion equation.  _Phys. Rev. Lett._ 86, 2712 (2001).
  * (37) Min, W., Luo, G., Cherayil, B. J., Kou, S. & Xie, X. S.  Observation of a power-law memory kernel for fluctuations within a single protein molecule.  _Phys. Rev. Lett._ 94, 198302 (2005).
  * (38) Kou, S. & Xie, X. S.  Generalized langevin equation with fractional gaussian noise: subdiffusion within a single protein molecule.  _Phys. Rev. Lett._ 93, 180603 (2004).
  * (39) Wei, Q.-H., Bechinger, C. & Leiderer, P.  Single-file diffusion of colloids in one-dimensional channels.  _Science_ 287, 625–627 (2000).
  * (40) Turiv, T. _et al._ Effect of collective molecular reorientations on brownian motion of colloids in nematic liquid crystal.  _Science_ 342, 1351–1354 (2013).
  * (41) Ochab-Marcinek, A. & Hołyst, R.  Scale-dependent diffusion of spheres in solutions of flexible and rigid polymers: mean square displacement and autocorrelation function for fcs and dls measurements.  _Soft Matt._ 7, 7366–7374 (2011).
  * (42) Cutland, N. J., Kopp, P. E. & Willinger, W.  Stock price returns and the joseph effect: a fractional version of the black-scholes model.  In _Seminar on stochastic analysis, random fields and applications_ , 327–351 (Springer, 1995).
  * (43) Burnecki, K. _et al._ Universal algorithm for identification of fractional brownian motion. a case of telomere subdiffusion.  _Biophys. J._ 103, 1839–1847 (2012).
  * (44) Ernst, D., Hellmann, M., Köhler, J. & Weiss, M.  Fractional brownian motion in crowded fluids.  _Soft Matt._ 8, 4886–4889 (2012).
  * (45) Weiss, M.  Single-particle tracking data reveal anticorrelated fractional brownian motion in crowded fluids.  _Phys. Rev. E_ 88, 010101 (2013).
  * (46) Mason, T. G. & Weitz, D.  Optical measurements of frequency-dependent linear viscoelastic moduli of complex fluids.  _Phys. Rev. Lett._ 74, 1250 (1995).
  * (47) Levernier, N., Bénichou, O., Voituriez, R. & Guérin, T.  Kinetics of rare events for non-markovian stationary processes and application to polymer dynamics.  _Phys. Rev. Res._ 2, 012057 (2020).
  * (48) Yaglom, A. M.  Correlation theory of processes with random stationary n th increments.  _Matematicheskii Sbornik_ 79, 141–196 (1955).
  * (49) Gripenberg, G. & Norros, I.  On the prediction of fractional brownian motion.  _J. Appl. Prob_ 33, 400–410 (1996).
  * (50) Anh, V. V. & Inoue, A.  Prediction of fractional brownian motion with hurst index less than 1/2.  _Bulletin of the Australian Mathematical Society_ 70, 321–328 (2004).
  * (51) Inoue, A. & Anh, V. V.  Prediction of fractional processes with long-range dependence.  _Hokkaido Mathematical Journal_ 41, 157–183 (2012).
  * (52) Eaton, M. L.  _Multivariate Statistics, A Vector Space Approach_ , vol. 53 (Institute of Mathematical Statistics Beachwood, Ohio, USA, 1983).
  * (53) Krapivsky, P., Mallick, K. & Sadhu, T.  Dynamical properties of single-file diffusion.  _J. Stat. Mech.: Theor. Exp._ 2015, P09007 (2015).
  * (54) Blumen, A., Von Ferber, C., Jurjiu, A. & Koslowski, T.  Generalized vicsek fractals: Regular hyperbranched polymers.  _Macromolecules_ 37, 638–650 (2004).
  * (55) Hosking, J. R.  Modeling persistence in hydrological time series using fractional differencing.  _Water resources research_ 20, 1898–1908 (1984).
  * (56) Dieker, T.  _Simulation of fractional Brownian motion_.  Master’s thesis (2004).
  * (57) Sokolov, I. M.  Cyclization of a polymer: first-passage problem for a non-markovian process.  _Phys Rev Lett_ 90, 080601 (2003).