Interlacing Relaxation and First-Passage Phenomena in Reversible Discrete and Continuous Space Markovian Dynamics
David Hartich Aljaz Godec

 Interlacing Relaxation and First-Passage Phenomena]Interlacing Relaxation and First-Passage Phenomena in Reversible Discrete and Continuous Space Markovian Dynamics


 
Mathematical Biophysics Group, Max-Planck-Institute for
  Biophysical Chemistry, Göttingen 37077, Germany
mailto:david.hartich@mpibpc.mpg.dedavid.hartich@mpibpc.mpg.de,mailto:agodec@mpibpc.mpg.deagodec@mpibpc.mpg.de


We uncover a duality between relaxation and first passage processes in ergodic reversible Markovian dynamics in both discrete and continuous state-space. The duality exists in the form of a spectral interlacing – the respective time scales of relaxation and first passage are shown to interlace. Our canonical theory allows for the first time to determine the full first passage time distribution analytically from the simpler relaxation eigenspectrum. The duality is derived and proven rigorously for both discrete state Markov processes in arbitrary dimension and effectively one-dimensional diffusion processes, whereas we also discuss extensions to more complex scenarios. We apply our theory to a simple discrete-state protein folding model and to the Ornstein-Uhlenbeck process, for which we obtain the exact first passage time distribution analytically in terms of a Newton series of determinants of almost triangular matrices.









§ INTRODUCTION



In his seminal work  <cit.> Kramers analyzed the kinetics of
chemical reactions in terms of diffusive barrier crossing, assuming that
the kinetic rate of a chemical reaction
corresponds to the inverse of the mean first crossing time. 
Ever since, first passage theory is at the heart of theoretical
descriptions of kinetics of chemical reactions 
<cit.>; see e.g. <cit.> for comprehensive reviews.
 
In a broader context, first passage concepts were invoked
in studies of kinetics in complex media, such as
reactions in fractal-like <cit.> and planar domains <cit.>,
in inhomogeneous cellular environments
<cit.>, in the study of
neural networks <cit.>, ultra cold atoms <cit.>,
 as well as in diverse narrow escape problems
<cit.> and so-called intermittent search strategies involving searching agents with internal dynamics
<cit.> (see also <cit.> for a review). 


First passage times play an important role in quantifying persistence properties
in non-equilibrium interacting many-body systems
<cit.>. More recent applications of first passage
concepts also include stochastic thermodynamics <cit.>,
in particular, fluctuation relations for stopping time statistics and stochastic entropy production in
driven molecular systems
<cit.> and in stochastic resetting processes <cit.>,
as well as uncertainty relations for first passage time statistics of fluctuating currents
<cit.> (see also <cit.>).

Moreover, our current understanding of the speed and precision of
transcription regulation in biological cells, and in particular of the
role of the so-called proximity effect in the co-regulation of genes,
<cit.> builds on first passage time ideas. 
The corresponding physical principles underlying these proximity effects
were explained in <cit.>. Notably, these works
revealed the inherent insufficiency of the mean first passage time
and traditional rate-based concepts  for a quantitative description of biophysical dynamics
in the so-called few encounter limit <cit.>. As a result,
a quantitative understanding of phenomena such as gene regulation <cit.> and the misfolding-triggered
pathological aggregation of proteins  <cit.>, which
are discussed in more detail in a related study <cit.>, requires the consideration of
the full statistics of first passage time.

Existing studies of the full first passage statistics in physical systems
typically focus on systems with continuous state-space dynamics,
whereas much less emphasis is put on discrete-space dynamics
<cit.>. Recent investigations of such discrete-state dynamics
include, for example, simple models of enzyme kinetics <cit.>
and novel numerical approximation schemes for studying  first-passage statistics based on Bayesian inference
 <cit.> (see also <cit.> for a recent review).  

 
Complementary to first passage processes are relaxation dynamics,
which by contrast do not terminate upon reaching a given threshold
for the first time. Relaxation phenomena in reversible diffusive
dynamical systems are nowadays well understood in terms of the eigenmodes and eigenvalues
of the underlying Fokker-Planck operators, which provide a generic and
very intuitive understanding of the dynamics of complex stochastic systems
<cit.>. Conversely, despite for allowing an
analogous spectral representation, a similar intuitive understanding of the
full first passage statistics and its physical implications remains
elusive. Notwithstanding, an important approximate link between 
the mean first passage time for escaping the deepest potential
basin and the corresponding slowest relaxation mode in the potential
was established in the seminal works of Matkovsky and Schuss
<cit.>, which has ever since been used routinely in
explaining relaxation phenomena in condensed matter
systems. Nevertheless, a deeper and more generic connection between
the two paradigms to date was not yet established. 

Here, we present the complete duality between relaxation and first
passage phenomena, which holds for all ergodic Markov
processes obeying detailed balance in both, continuous and discrete 
state-space, in which the absorbing target is effectively one-dimensional. The duality emerges in the form of a spectral
interlacing, which we prove rigorously by combining spectral-theoretic, matrix-algebraic
and Greens function-theoretic concepts.
On the one hand the duality
allows for an intuitive generic understanding of first passage
phenomena in terms of relaxation eigenmodes. On the other hand, it
enables us to determine the full first passage time statistics exactly
from the corresponding relaxation eigensystem.  The formalism is exact
and holds for all reversible Markovian systems governed by a master
equation in arbitrary dimensions or by a Fokker-Planck equation, and
therefore unifies the theoretical treatment of discrete and
continuous space phenomena. 
We note the spectral interlacing in the case of a discrete state dynamics
has also recently be deduced from a `lumping' of the state dynamics <cit.>.

To illustrate the predictive power of the formalism in practice, we
here predominantly focus on systems with discrete state-space dynamics,
whereas continuous space dynamics are treated in more detail in a related study <cit.>. In particular, we here apply our theory to a simple discrete-state protein folding model
and to diffusion in a harmonic potential, also know as the
Ornstein-Uhlenbeck process. Notably, we obtain, to the best of our
knowledge, for the first time an exact analytical solution for the full first passage time distribution of the Ornstein-Uhlenbeck process in the time domain.

The paper is organized as follows. In Sec. <ref>
we present a canonical formulation the first passage problem
applicable to both discrete states-pace and continuous Fokker-Planck
dynamics. Sections <ref> and <ref>
provide a step-by-step explanation of how one can exactly determine
the first-passage distribution from the corresponding relaxation process, and also contain rigorous proofs of the duality in
discrete and continuous state-space dynamics, respectively.
We apply the duality framework in Sec. <ref> to determine the
first passage statistics for a simple
protein folding model and for the Ornstein-Uhlenbeck process. A concluding
perspective is provided in Sec. <ref>.
In <ref> we derive a compact representation of the long-time
asymptotics of the first passage time distribution, which inter alia
extends our results for the long time asymptotics from equilibrium systems to irreversibly driven systems.




§ FUNDAMENTALS




 §.§ Relaxation and first passage

We assume that the probability density to find the system in state
x at time t upon evolving from an initial state x_0 according to
microscopically reversible dynamics, P(x,t|x_0), is governed
by

    _t P(x,t|x_0)= P(x,t|x_0),

where  is a linear reversible operator, which will be specified
below. We consider two classes of operators: (DS) discrete state Markov jump
process, where x and x_0 assume only a finite number of states,
and (FP) continuous Markovian diffusion governed by a Fokker-Planck equation.


For discrete Markov state models of class (DS)
the dynamics is governed by

    P(x,t|x_0)≡∑_x'=0^ L(x,x') P(x',t|x_0),

where x,x'=0,1,…, denote the discrete states, L(x,x') is
the rate of jumping from state x' to state x  (x≠ x') and
-L(x,x)=∑_x'≠ x L(x',x) is the total rate of leaving state x guaranteeing conservation
of probability (∑_x_t P(x,t|x_0)=0).
In order to have reversible dynamics
we need to
additionally impose detailed balance, i.e. the constraint
L(x,x')/L(x',x)=exp[β U(x')-β U(x)] (see., e.g.,
<cit.>), which assures that the system will relax to a
Boltzmann distribution in a potential U(x) on ergodic timescales
P_eq(x)∝^-β U(x), where β=1/k_B T is the
inverse thermal energy.
We call such a reversible ergodic process that conserves probability a relaxation process.
If we add an absorbing point at x=a we call the resulting process
a first passage process or in short absorption, which we introduce in the following way.
First, we modify the generator (→) such that all transitions corresponding
to jumps out of the absorbing state a are removed, i.e., the
elements of the
first passage generator read

    L_a(x,x')={[          0   if x'=a,;    L(x,x') otherwise. ].

Using a bra-ket matrix notation <cit.> we rewrite this equation as

    =-|a⟩⟨a|,

where |a⟩=(0,… ,0,1,0,…,0)^=(|a⟩)^ is a vector with all entries except the ath one; consequently,
we identify L_a(x,x')=⟨x||x'⟩.
The first passage time density to reach state a at time t
starting from x_0
is then formally defined by

    ℘_a(t|x_0)=_t⟨a|^ t|x_0⟩=⟨a|^ t|x_0⟩,

which is nothing but the normalized probability flux into state a
with
∫_0^∞℘_a(t|x_0) t=1.
Note that with Eq. (<ref>) we use the convention that
|a⟩ is the unique stationary solution with |a⟩=0.


For a continuous space Markovian diffusion the transition probability density
function (the `propagator') instead obeys the Fokker-Planck equation (<ref>) 

    P(x,t|x_0)   =-_x j(x,t|x_0)
       ≡_x D[β U'(x)+_x]P(x,t|x_0),

where j(x,t|x_0) is the probability current, D is the diffusion
constant, -U'(x)=-_xU(x) is a force field generated by the potential U(x)
at position x, and β is the inverse temperature, which we set to β≡ 1 to express energies in units of T from now on.
The scenario with reflecting barriers at x=b_± with j(b_±,t|x_0)=0
<cit.>, we term a relaxation process, where
b_±=±∞ correspond to so-called natural boundary conditions [For natural boundary conditions the current and the probability density both vanish, i.e., lim_b→±∞j(b,t|x_0)=lim_b→±∞P(b,t|x_0)=0.
].

Conversely, an absorbing boundary at x=a enters the Fokker-Planck equation via the Dirichlet boundary condition P(a,t|x_0)=0,
without altering the partial differential equation (<ref>), i.e.,
the first passage operator still reads =_xD[U'(x)+_x].
However, here the first passage time density becomes the probability flux into state a.
For convenience we use the operator  as shorthand for
Eq. (<ref>) under the boundary condition
P(a,t|x_0)=0.



We note that without an absorbing point both, dynamics governed by the master equation (<ref>) and the Fokker-Planck equation (<ref>)
relax to the Boltzmann distribution P_eq(x)∝exp[-U(x)], whereas
with the absorbing boundary condition the particle will eventually reach the target
with probability 1.





 §.§ Eigendecomposition


Since  is assumed to generate a reversible Markov process,
we can expand the generator
 in a bi-orthogonal eigenbasis <cit.>. Denoting the
eigenvalues of the relaxation process by λ_k and the
corresponding left (right) eigenvectors by ⟨_k|
(|_k⟩), respectively, the generators from
Eqs. (<ref>) and (<ref>) 
become in the respective eigenbases

    =-∑_k λ_k |_k⟩⟨_k|,

where λ_0=0≤λ_1≤…, and ⟨_k|_l=δ_kl.
We assume the eigenvalues to be ordered
such that λ_k≤λ_k+1, and the generator to be irreducible
λ_0=0<λ_1, which means that there is a unique equilibrium state <cit.>.
Note that for a Fokker-Planck
equation with reflecting barriers at x=b_± (relaxation) the eigenfunction _k(x)≡⟨x|_k must satisfy the zero flux condition -D[U'(x)+_x]_k(x)|_x=b_±=0, with ⟨x|_0∝^-U(x).

The generator with the absorbing point at state a, can similarly be expanded in a bi-orthogonal set of eigenfunctions

    =-∑_k μ_k |_k⟩⟨_k|,

where μ_k is the k-th eigenvalue and ⟨_k| (|_k⟩) denote the corresponding left (right) eigenfunctions
of the first passage process.
Without loss generality we use an ordered labeling such that μ_k≤μ_k+1,
where 0<μ_1.





The left and right eigenvectors of the absorption (at position x≠
a) as well as of the relaxation process
are related via ⟨x|_k∝^U(x)⟨x|_k and
⟨x|_k∝^U(x)⟨x|_k, respectively.
In the case of a discrete number of states, the lowest eigenvalue of the generator (<ref>)
will be μ_0=0 with the right eigenfunction |_0⟩=|a⟩,
whereas for Fokker-Planck dynamics one imposes the boundary condition
⟨a|_k=0.

In a previous work an explicit Newton series expression for μ_1 in
terms of a series of almost triangular matrices was derived <cit.>, which corresponds to a large deviation limit t→∞.
One of our main goals here is to obtain the full first passage statistics ℘_a(t|x_0) explicitly
in terms of relaxation eigenmodes. Our theory builds on the renewal
theorem, which we briefly review
in the following subsection.






 §.§ Renewal theorem

The classical renewal theorem provides
a well known implicit connection between first passage
and relaxation processes. It relates the 
probability density of the freely propagating system to be in state x at time t
upon starting from a state x_0, to the first passage distribution
℘_a(t|x_0) from x_0 to a:

    P(x,t|x_0)=∫_0^tτ  P(x,t-τ|a) ℘_a(τ|x_0),

where both P(x,t|x_0) and ℘_a(t|x_0) admit a spectral representation

    P(x,t|x_0)=⟨x|^ t|x_0⟩=∑_k
     ⟨x|_k⟨_k|x_0^-λ_kt

and

    ℘_a(t|x_0)=∑_k≥ 1 w_k(x_0)μ_k^-μ_kt ,

respectively. In other words, a system starting from state x_0 must pass through state a
before reaching  the final state x,
which for an effectively 1-dimensional Fokker-Planck necessarily means x_0< a≤ x or x_0> a≥ x. In
Eq. (<ref>) we introduced in the first passage weights

    w_k(x_0)=
    {[                  -⟨a|_k⟨_k|x_0                        for DS,; -σ_± D⟨x|_k/ x|_x=a⟨_k|x_0/μ_k                        for FP, ].

for discrete state (DS) and Fokker-Planck (FP) dynamics, respectively, which
must satisfy ∑_kw_k=1 with the first nonzero weight being
strictly positive w_1(x_0)>0, and where we introduced σ_±≡sign(a-x_0).
Note that the first line of Eq. (<ref>), i.e. the DS case,
is equivalent to w_k(x_0)=∑_x≠ a⟨x|_k⟨_k|x_0
with _a from Eq. (<ref>). 
In the case of FP dynamics the second line of
Eq. (<ref>) is equivalent to
 w_k(x_0)≡∫⟨x|_k⟨_k|x_0 x,
which follows from a partial integration using both
Eq. (<ref>) and Eq. (<ref>).



In the case of x=a the renewal theorem (<ref>) has the simple interpretation:
a system being in state a at time t must have arrived at that point
at some earlier time τ for the first time  (τ≤ t),
and then returned to the same position again at time t, where τ=
t corresponds to the time of first arrival.

Laplace transforming the renewal theorem (<ref>), where a generic function f is transformed according to f̃(s)≡∫^-s tf(t) t, 
yields <cit.>

    ℘̃_a(s|x_0)=P̃(x,s|x_0)/P̃(x,s|a)
     =∑_k(s+λ_k)^-1⟨x|_k⟨_k|x_0/∑_k(s+λ_k)^-1⟨x|_k⟨_k|a.

Based on this well known renewal theorem we construct in the following
section a method that allows to determine explicitly the first passage time statistics ℘_a (t|x_0) exactly in terms of the relaxation process, i.e., we
render Eq. (<ref>) explicit in the time domain.



§ PRINCIPAL RESULT FOR DISCRETE STATE SYSTEMS


Starting from the renewal theorem (<ref>), we now derive an expression for the first passage time density
for discrete state Markov processes
in terms of
relaxation modes in the following three steps.
The first step involves a crucial relation
between the eigenvalues of the
relaxation process λ_k and absorption process μ_k, which
are here shown to interlace

    λ_k-1≤μ_k≤λ_k

for k=1,…,.
For effectively one dimensional finite lattice models
with the target at an outer edge
these inequalities become strict

    λ_k-1<μ_k<λ_k,

which will also apply identically to Fokker-Planck dynamics discussed in Sec. <ref> in which case we formally assume =∞.
In the second step we exactly express the first
passage eigenvalues μ_k in the form of a Newton series of
determinants of
almost triangular matrices, which generalizes the result for the slowest mode μ_1 from <cit.> to all first passage modes.
The third and final step corresponds to a straightforward application
of the residue theorem, which is used to determine the first passage
weights w_k(x_0).





 §.§ Interlacing of eigenmodes (step 1)



For a discrete system with +1 states the eigenvalues λ_k and μ_k
correspond to the roots of the respective characteristic polynomials

    χ(s)   ≡(𝕀 s-)=s∏_i=1^(s+λ_i),
    χ_a(s)   ≡(𝕀 s-)=s∏_i=1^(s+μ_i),

i.e., χ(-λ_k)=0 and χ_a(-μ_k)=0.
Inserting Eq. (<ref>), which is =-|a⟩⟨a|, into the
second characteristic polynomial (<ref>)
and using the matrix
determinant lemma establishes a link between the two characteristic polynomials

    χ_a(s)=χ(s)+⟨a|(𝕀 s-)|a⟩,

where (A) is called the adjugate of a matrix A satisfying
Cramer's rule A(A)=(A)𝕀. We note that
the same mathematical concepts have been used  recently to determine
the stalling distribution of irreversibly driven systems
(cf. deletion-contraction formula in <cit.>).

The adjugate of a diagonal matrix D
with elements D_ii=d_i (D_ij=0 if i≠ j) is diagonal as well, with elements (D)_ii=∏_j≠ id_j. Consequently, the bi-orthogonal
expansion (<ref>)
implies

    (𝕀 s-)=∑_i=0^|_i⟩⟨_i|∏_j=0
    j≠ i^(s+λ_j),

which inserted into Eq. (<ref>)
gives

    χ_a(s)   =χ(s)-∑_i=0^⟨a|_i⟨_i|aλ_i∏_j=0
    j≠ i^(s+λ_j)
     ,

where we used the eigenvalue equation
⟨_i||a⟩=-λ_i⟨_i|a. Eq. (<ref>)
constitutes an essential step in our calculations, which allows us to
express the diagonal of the relaxation propagator P̃(a,s|a)
solely in terms of eigenvalues μ_k and λ_k (see the following subsection for more details).

Moreover,
the characteristic polynomials of the first passage process χ_a(s) and relaxation process χ(s) change
sign one after the other, since detailed balance imposes ⟨a|_i⟨_i|a≥0 for all i=0,…,,
which proves that the eigenvalues of the first passage process μ_k and eigenvalues of the relaxation process λ_k
interlace according to Eq. (<ref>).
We note that this result is directly related to the interlacing of eigenvalues generated from a  “lumping” of
states which is proven in <cit.>.
In the following paragraph we briefly discuss the scenario, in which
the interlacing of eigenvalues becomes strict (<ref>),
which will be the case for systems with
Fokker-Planck dynamics discussed in Sec. <ref>.





The stronger condition  (<ref>) holds if all eigenfunctions are nonzero at the target |⟨a|_k|>0 and all relaxation eigenvalues are non-degenerate, that is, λ_i-1<λ_i for all i=1,…,.
One can show that this condition is always trivially satisfied for 1-dimensional models (L(x,x')=0 if |x-x'|>1), in which the target a is placed at the border (e.g., a= or a=0); see inset
of Fig. fig:illustration_proofa for such an exemplary 3-state system.

Inserting the relaxation eigenvalues s=-λ_k into the characteristic polynomial
of the first passage process (<ref>)
yields

    χ_a(-λ_k)
        =(-1)^k+1λ_k
     ⟨a|_i⟨_i|a∏_i=0
    i≠ k^ |λ_i-λ_k|,

where we used the relations
λ_i-λ_k< 0 for all i<k and
λ_i-λ_k>0 for all i>k, as well as χ(-λ_k)=0.
Since for k≥ 1 each eigenvalue is positive (λ_k>0),
the characteristic polynomial of the first passage process χ_a(-λ_k)
is
equal to (-1)^k-1 multiplied by a positive constant.
Consequently, 
χ_a changes sign exactly once between any two consecutive relaxation modes -λ_k<s<-λ_k-1.
The fact that χ_a and χ are polynomials of the same degree +1 forbids
more than a single root,
and hence implies the strict interlacing of eigenvalues from
Eq. (<ref>), which completes the proof.
The aforementioned reasoning is illustrated in Fig. fig:illustration_proofa for a simple three state model
in which the vertical arrows represent Eq. (<ref>).

For fine-tuned systems in which the target is not located at the very outer position
(see e.g., Fig. fig:illustration_proofb) or systems that are not effectively one dimensional
(see e.g., Fig. fig:illustration_proofc) the strict interlacing theorem (<ref>) can be violated,
whereas the “slightly weaker” interlacing condition (<ref>)  still holds.





 §.§ Diagonal of the relaxation propagator in terms of bare eigenvalues


Using the results from the previous subsection we are now in the position to
represent P̃(a,s|a) (i.e. Eq. (<ref>) with x_0=x=a),
using only the eigenvalues of both the first passage process and the
relaxation process, μ_k and λ_k, respectively.
Laplace transforming the eigenmode expansion in Eq. (<ref>) assuming x_0=x=a yields

    P̃(a,s|a)
     =P_eq(a)/s+∑_k=1^⟨a|_k⟨_k|a/s+λ_k,

where we identified the equilibrium probability density ⟨a|_0⟨_0|a=P_eq(a) in the first term.
Comparing now P̃(a,s|a) in (<ref>)
with χ_a from Eq. (<ref>) and χ from Eq. (<ref>)
yields after some algebra

    P̃(a,s|a)=χ_a(s)/sχ(s)=1/s∏_i=1^(s+μ_i)/(s+λ_i).

The second equality in Eq. (<ref>) follows from Eq. (<ref>).
Hence, P̃(a,s|a) encodes
the eigenvalues of both, the relaxation and the first passage processes.
Due to Eq. (<ref>) P̃(a,s|a) contains only simple poles
and decreases monotonically in s between any two consecutive poles
since ⟨a|_k⟨_k|a≥0.
If ⟨a|_k⟨_k|a>0 (e.g., 1d models
with the target at at the border) each root λ_k of P̃(a,s|a) represents a first passage eigenvalue s=-μ_k (k≥ 1), which is located in between two relaxation modes
λ_k-1<μ_k<λ_k, thus providing an alternative proof of relation (<ref>) <cit.>.
In the following section we determine the roots of the diagonal of the propagator explicitly, which due to
Eq. (<ref>) correspond to first passage eigenvalues μ_k.

Let us briefly reformulate P̃(a,s|a) 
in a way that  can also be applied to continuous systems with an infinite number of states.
Isolating the equilibrium probability, which is the first term in Eq. (<ref>), from the product formula (<ref>)
yields

    P̃(a,s|a)=P_eq(a)/s∏_k=1^(1+s/μ_k)/(1+s/λ_k).

Since μ_k,λ_k increase monotonically with k we will later be able to
adopt these results to systems governed by Fokker-Planck dynamics, which formally
corresponds to the limit →∞ for which the product in
Eq. (<ref>) still
converges.




 §.§ From the relaxation spectrum to the first passage time spectrum (steps 2 and 3)


Based on the interlacing theorem presented in Eq. (<ref>), which is
also given in Eq. (4) in a related work <cit.>, we can determine the full first passage time spectrum
{μ_k,w_k(x_0)} from the corresponding relaxation spectrum,
{λ_k,|_k⟩,⟨_k|}.
For simplicity we first consider the eigenvalues to be both ordered
λ_k<λ_k+1 and non-degenerate, and also assume that ⟨a|_k⟨_k|a>0 holds for all values of k.
The extension to situations with ⟨a|_k⟨_k|a=0, which also includes degenerate eigenvalues, for some k is straightforward and will be dealt with at the end of this subsection.


Before determining the weights w_k, we first determine the first passage
eigenvalues μ_k, which were shown to be encoded in the roots of P̃(a,s|a) in Eq. (<ref>).
We introduce the k^*th “modified diagonal of the propagator”

    F_k^*(s)   ≡ (s+λ_k^*)P̃(a,s|a)
       =⟨a|_k^*⟨_k^*|a+(s+λ_k^*)∑_l=0
    l≠ k*^⟨a|_l⟨_l|a/s+λ_i,

which still encodes all of the first passage eigenvalues {μ_k}
according to Eq. (<ref>), i.e., it has exactly the same roots
as P̃(a,s|a). However, in contrast to P̃(a,s|a) the modified function F_k^*(s) is strictly concave within the interval -λ_k^*+1< s<-λ_k^*-1,
which can easily be confirmed by taking the second derivative and realizing that F̈_k^*(s)≡_s^2F_k^*(s)<0 holds within the region of interest -λ_k^*+1< s<-λ_k^*-1.

For k^*=k and k^*=k-1 the modified functions   F_k(s) and F_k-1(s) both are strictly concave within the interval -λ_k<s<-λ_k-1
and, consequently, also locally concave around the kth first passage eigenvalue s=-μ_k, i.e., F̈_k(-μ_k) and F̈_k-1(-μ_k)<0.
Moreover, both functions F_k(s) and F_k-1(s) allow a Taylor
expansion around the  midpoint μ̅_k≡(λ_k+λ_k-1)/2
that converges within the whole interval -λ_k<s<-λ_k-1
including the root s=-μ_k at which F_k(-μ_k)=F_k-1(-μ_k)=0.

The method we present in the following is an analytical technique based on the principles of Newton iteration, which is a
simple root finding algorithm that is guaranteed to work for functions that are both negative and concave between the starting point and the first root.
Hence, to determine the kth eigenvalue we accordingly choose the
modified function 

    f(s,k)= F_k^*(s),

such that

    k^*=
     {[                k if F_k(-μ̅_k)<0,;              k-1       otherwise. ].
     ,

which guarantees both negativity f(s,k)≤ 0 and concavity _s^2f(s,k)≤0 between
s=-μ̅_k and s=-μ_k.

According to the interlacing theorem (<ref>) s=-μ_k is the only zero f(-μ_k,k)=0 within the interval -λ_k<s<-λ_k-1.
With the midpoint starting condition μ̅_k=(λ_k+λ_k-1)/2
the kth first passage eigenvalue can be represented exactly in
a series of determinants of almost triangular matrices

    μ_k=μ̅_k+∑_n=1^∞f_0(k)^n/f_1(k)^2n-1_n(k)/(n-1)!,

where f_n(k) is the nth derivative of f(s,k) as defined in (<ref>) with respect to s at s=-μ̅_k, and
_n(k) stands for an almost triangular matrix 
with elements <cit.>

    _n^i,j(k)=f_i-j+2(k)Θ(i-j+1)/(i-j+2)![n(i-j+1)Θ(j-2)
     +iΘ(1-j)+j-1],

with Θ(l) denoting the Heaviside step function (Θ(l)=1 if l≥0)
and i,j=1,2,…,n-1. Moreover, we adopt the convention _1(k)=1.
We note that this method generalizes 
the method recently derived to determine the slowest first passage mode μ_1 <cit.> to all first passage eigenmodes μ_k.

Let us briefly repeat the two crucial steps towards Eq. (<ref>). First, the interlacing theorem (<ref>) guarantees that the Taylor series 
f(s,k)=∑_i f_i(k)(s+μ̅_k)^i around the midpoint
μ̅_k=(λ_k+λ_k-1)/2 converges in the entire spectral interval 
-λ_k<s<-λ_k-1,
which also contains the first passage eigenvalue s=-μ_k. Second,
due to F_k^*(s) in
Eqs. (<ref>)-(<ref>)
the function f(s,k) is strictly concave and negative between s=-μ̅_k and s=-μ_k, which in turn guarantees the convergence of the
explicit Newton series (<ref>).

Eqs. (<ref>)-(<ref>) provide a universal method for determining
explicitly first passage eigenvalues from the corresponding relaxation spectrum
and constitute the central result of this work. We show in the <ref> a simpler derivation of μ_1 as well as a compact approximation of the principal first passage eigenvalue μ_1, which is particularly useful in the case of time scale separation μ_1≪λ_1 (or λ_1≪λ_2).
Furthermore, <ref> provides a generalization of the long time
asymptotics from systems with reversible dynamics to irreversibly
driven systems.

In the following we briefly comment on the practical
implementation of the exact result for μ_k  to render Eqs. (<ref>)-(<ref>) fully explicit. The weights w_k will be determined afterwards in this subsection.
The nth derivative of F_k^*(s) with respect to s
at s=-μ̅_k,  f_n(k)≡∂_s^nf(s,k)|_s=-μ̅_k,
can be written explicitly as

    f_0(k)=⟨a|_k^*⟨_k^*|a+∑_ l|l≠ k^*⟨a|_l⟨_l|a(μ̅_k-λ_k^*)/(μ̅_k-λ_l),
       f_n≥ 1(k)=n!∑_ l|l≠ k^*⟨a|_l⟨_l|a(λ_l-λ_k^*)/(μ̅_k-λ_l)^n+1,

where k^*=k or k^*=k-1 is chosen according to Eq. (<ref>). Note that condition (<ref>) is equivalent
to the condition f_0(k)≤ 0, implying the first line of Eq. (<ref>)
to be either negative for k^*=k or for k^*=k-1,
i.e. one has to evaluate the first line of Eq. (<ref>) for k^*=k:
if f_0(k)>0 one must to change k^* to k^*=k-1 and reevaluate f_0(k).
Once one has determined k^*(k) and f_0(k) one can proceed with the second line of Eq. (<ref>) to determine f_n(k)/n! and insert the result in the almost triangular matrix (<ref>).
The determinant of almost triangular matrices can be
calculated elegantly using the simple recursion relation from <cit.>, see also <cit.> for an efficient numerical implementation.



Having obtained the first passage eigenvalues,
the weights of the first passage time distribution can be calculated
using the standard residue theorem. The Laplace transform of the spectral expansion of the first passage time density (<ref>) reads

    ℘̃_a(s|x_0)≡∑_k w_k(x_0)μ_k/s+μ_k.

Using the residue theorem to invert the Laplace transformed renewal theorem
(<ref>)
yields

    w_k(x_0)   =P̃(a,-μ_k|x_0)/μ_kṖ̃̇(a,-μ_k|a)
       =∑_l(1-λ_l/μ_k)^-1⟨a|_l⟨_l|x_0/∑_l(1-λ_l/μ_k)^-2⟨a|_l⟨_l|a,

where Ṗ̃̇(a,s|a)=_sP̃(a,s|a) is taken at s=-μ_k. The explicit Newton series (<ref>) along with the first passage weights (<ref>) fully characterize the first passage time
distribution ℘_a(t|x_0)=∑_k w_k(x_0)μ_k^-μ_k t in terms of relaxation
eigenmodes {λ_k,_k}. This completes our third and final step, which allows, for the first time, to analytically deduce first passage time statistics directly from relaxation eigenmodes. We call this relation the explicit forward duality between first passage and relaxation. This completes the central result of this paper.


The spectral representation is very useful for determining the moments
of the first passage time, ⟨ t^n⟩≡∫
t^n℘_a(t|x_0) t=n!∑_kw_k(x_0)μ_k^-n.
Moreover, as explained in more detail in a related work <cit.>, the
full spectral expansion is required for a correct explanation of
kinetics in the 
so-called few encounter limit, where N molecules
starting from position x_0 are searching for the target at a. 
The probability density that the first molecule out of N arrives at
time t at a for the first time for this case
becomes ℘_a^(N)(t|x_0)=N℘_a(t|x_0)[∫_t^∞℘_a(τ|x_0)]^N-1τ, which can be understood as follows.
The probability that the first N-1 molecule have not yet reached the target will be given by [∫_t^∞℘_a(τ|x_0)τ]^N-1, while the Nth particle
arrives at a with a rate ℘_a(t|x_0); hence the probability density that any particle
out of N molecules arrives at the target for the first time
according to ℘_a^(N)(t|x_0). Further details of the N-particle
problem and in particular the physical implications of the
few-encounter limit are discussed in a related study
<cit.>.


Let us now briefly generalize the method to systems with degenerate eigenvalues
or vanishing relaxation modes.
An eigenfunction that vanishes at the target ⟨a|_k=0
will have a vanishing spectral weight as a result of
Eq. (<ref>). Hence, `manually' removing such modes will not
affect the first passage time distribution ℘_a. Moreover, if a relaxation eigenvalue λ_k 
is degenerate we define 

    Ψ_k(a,x_0)≡∑_k'|λ_k'=λ_k⟨a|_k'⟨_k'|x_0,

and replace ⟨a|_k⟨_k|x_0→Ψ_k(a,x_0)
as well as ⟨a|_k⟨_k|a→Ψ_k(a,a)
and take the sums in Eq. (<ref>) over all different values of λ_k.
After renumbering all distinct contributing eigenvalues we obtain a strict
interlacing (<ref>). Therefore, we can apply our standard forward duality also to degenerate eigensystems.
In the next subsection we will briefly derive a formal backward duality
after which we reformulate the results from this subsection to
continuous Fokker-Planck dynamics.



 §.§ Backward duality


In contrast to the explicit forward duality, which was presented in the previous subsection, an explicit reverse relation in the time-domain could not be established.
In Laplace space, however, the forward duality can be inverted to give
a backward duality as follows. Inserting the first passage generator
=-|a⟩⟨a| from (<ref>) into the Laplace transform
of the propagator P̃(a,s|x_0)=⟨a|(𝕀 s-)^-1|x_0⟩ and
using the Sherman-Morrison-Woodbury formula
yields

    P̃(a,s|x_0)=⟨a|(𝕀 s-)^-1|x_0⟩/1-⟨a|(𝕀 s-)^-1|a⟩.

Let us now insert the expression for the first passage time distribution from Eq. (<ref>), which can be written as ℘_a(s|x_0)=s⟨a|(𝕀 s-)^-1|x_0⟩, 
into Eq. (<ref>) to obtain

    P̃(a,s|x_0)=℘̃_a(s|x_0)/s/1 -∑_x=0^ L(x,a)℘̃_a(s|x)/s,

where L(x,a)=⟨x||a⟩ is the generator of the relaxation process.
Notably, this is expression corresponds to the backward duality
and is the formal inverse of the renewal theorem,
where ℘̃_a(s|x_0)/s is the Laplace transform of the cumulative first passage time distribution ∫_0^t℘_a(τ|x_0)τ.




§ PRINCIPAL RESULT FOR FOKKER-PLANCK DYNAMICS




 §.§ Greens function with natural boundaries

We restrict our discussion to effectively 1-dimensional dynamics, which include
diffusion in d dimensions in an isotropic potential as discussed in <cit.>, where d may also be fractal.
Introducing an absorbing target at position a splits the
first passage problem into two cases (I) x_0<a and (II) x_0>a.
Case (I) corresponds to an absorption from the left, and case (II) to an absorption from the right.
In the following paragraph we demonstrate that all first passage modes μ_k
of both distinct cases (I) and (II) are entirely encoded in  P̃(a,s|a), which allows to formulate the results from Sec. <ref>
also for systems with Fokker-Planck dynamics.


Laplace transforming the Fokker-Planck equation (<ref>) yields

    (-s)P̃(x,s|x_0)=-δ(x-x_0)

where =-_xD[β U'(x)+_x] and x_0 is the initial position of the relaxation process.
Eq. (<ref>) is a inhomogeneous linear differential equation
which can be solved using the standard Green's function approach.
First, we find the two independent solutions v_±(x,s) of the homogeneous problem (-s)v_±(x,s)=0, where we use the label “-” and “+”
for the solution satisfying the left and right boundary condition, respectively.
That is, a diffusion process within an interval b_-<x<b_+
imposes the probability current j_±(x,s)≡ - D[β U'(x)+_x]v_±(x,s)
to vanish at the boundaries, i.e. j_±(b_±,s)=0. The special case of so-called natural boundary conditions correspond to the limit  lim_x→±∞v_±(x,s)= 0
or analogously lim_x→±∞j_±(x,s)= 0, that is, b_±=±∞.
The full solution P̃(x,s|x_0) of (<ref>) is a continuous function in x with a discontinuity of its first derivative at x=x_0.
Using the scaled Wronskian
[For convenience we defined with the scaled Wronskian with the current function j_±(x,s) instead of the first derviative _xv_±(x,s), i.e., W_s(x)/D would represent the standard definition of the Wronskian.]

    W_s(x)
       ≡ D[v_-(x,s)_x v_+(x,s)-v_+(x,s)_x v_-(x,s)],
       =v_+(x,s)j_-(x,s)-v_-(x,s)j_+(x,s)
       = [ v_+(x,s) v_-(x,s); j_+(x,s) j_-(x,s) ]

the propagator, which satisfies the proper jump condition of the first derivative (current function) at x=x_0, becomes

    P̃(x,s|x_0)=
     {[ v_+(x,s)v_-(x_0,s)/W_s(x_0)                  if x_0≤ x,; v_-(x,s)v_+(x_0,s)/W_s(x_0)                  if x_0≥ x. ].

We note that the Wronskian (<ref>) is proportional to the Boltzmann factor (see, e.g., Ref. <cit.>), i.e., W_s(x)=W_s(x_0)exp[β U(x_0)-β U(x)]. Hence using the renewal theorem (<ref>)
and P̃(a,s|a) as well as P̃(a,s|x_0)
from Eq. (<ref>)
yields
the Laplace transform of the first passage time distribution

    ℘_a(s|x_0)=
     ^β U(x_0)-β U(a)×{[ v_-(x_0,s)/v_-(a,s)           if x_0<a,; v_+(x_0,s)/v_+(a,s)           if x_0>a. ].

The two independent functions v_± (x,s) are entire functions
without any poles in s <cit.>, and in turn encode in their
roots all first passage eigenvalues s=-μ_k. In particular v_- encodes all first passage modes
from case (I) x_0<a, and v_+ encodes all first passage modes from case (II), in which the particle is absorbed from the right x_0>a. Due to 
Eq. (<ref>) the zeros of P̃(a,s|a) at s=-μ_k determine the first passage spectrum.
Hence, all results from Sec. <ref> hold identically for continuous systems as well. However, the sums are here not finite, i.e., =∞. For example, P̃(a,s|a) becomes

    P̃(a,s|a)=∑_l=0^∞_l(a)_l(a)/s+λ_l,

where _l is the lth right eigenfunction of the Fokker-Planck operator
satisfying _l(x)=-λ_l_l(x) with the corresponding
left eigenfunction _l(a)∝^β U(a)_l(a) and normalization ∫_b_-^b_+_l(x)_l(x) x=1.
The first passage modes μ_k can then be determined with Eqs. (<ref>)-(<ref>), where k^* (k^*=k or k^*=k-1) must be chosen such that
f_0(k)<0  holds in Eq. (<ref>) with ⟨a|_l⟨_l|a≡_l(a)_l(a).
Concurrently, the first passage weights follow from Eq. (<ref>).

The formal backward duality from Sec. <ref>,
however, must be adopted as follows. After some tedious algebra we
obtain formally the exact inverse duality in the form of

    P̃(x,s|x_0)=σ_±^β U(x_0)-β U(x)℘_x_0(s|x)/D/ x_0ln[℘_x_0(s|x)℘_x(s|x_0)],

where sign σ_±=-1 if x_0< x and σ_±=+1 if x_0> x;
Eq. (<ref>)
can easily be verified by inserting
the  Wronskian (<ref>) and the first passage time distribution (<ref>) into the right hand side
of Eq. (<ref>), and comparing the result with the propagator from Eq. (<ref>). Notably,  this inverse duality is the continuous version of Eq. (<ref>).




 §.§ Relaxation under reflecting boundary conditions and strict
  spectral interlacing


In the previous subsection the target a divided the phase space into two regions,
which implies that the first passage modes for the cases (I) and (II) separate
into “left” and “right” modes as well. For example, if x_1<a and x_2>a
one of the first passage weights w_k(x_1) or w_k(x_2) must typically be zero for all values of k. 
If one uses just the first M modes to approximate the propagator [cf. Eq. (<ref>)]
the zeros of the right hand side of

    P̃^M(a,s|a)≡∑_k=0^M_k(a)_k(a)/s+λ_k,

become approximations of the first passage modes and, hence, the weights
w_k^M(x_1) and w_k^M(x_2) deduced from Eq. (<ref>) will only satisfy w_k^M(x_1)≪ w_k^M(x_2)
(or w_k^M(x_1)≫ w_k^M(x_2)) for finite M, i.e., modes from case (I) and (II)
mix. Such a mixing can be avoided entirely if the relaxation process is analyzed with
a reflecting boundary at the target position a (b_+=a or b_-=a).

The result for reflecting boundary conditions j_+(a)=0 (case (I)) and j_-(a)=0  (case (II)) 
is automatically obtained by the following replacement:

    v_±(x)   → v_±(x,s)j_∓(a,s)-j_±(a,s)v_∓(x,s),
    
      j_±(x)   → j_±(x,s)j_∓(a,s)-j_±(a,s)j_∓(x,s),

respectively, which inserted into the scaled Wronskian (<ref>)
at x=a
yields

    W_s(a)=
    {[     v_+(a,s)j_-(a,s)  if  case (I) x_0<a,;    -v_-(a,s)j_+(a,s) if  case (II) x_0>a. ].

Utilizing the Wronskian for the reflecting boundary condition in Eq. (<ref>) 
yields the diagonal of the propagator in the form

    P̃(a,s|a)≡lim_ϵ→0P̃(a±ϵ,s|a)=±v_±(a,s)/j_±(a,s)

For two linearly independent functions v_± with nonzero Wronskian (<ref>) the zeros of v_± and j_± are different. Hence, the zeros of v_- in s contain
only the first passage modes for the case (I) x_0<a, whereas the
zeros of v_- do not contain zeros of first
passage modes corresponding to the case (II) x_0>a.

Let us from now on just focus on a case (I), in which x_0<a, since case (II) follows by analogy.
For case (I) we consider the Fokker-Planck operator  from Eq. (<ref>)
with zero current condition at x=a and natural boundary condition for x→-∞.
To that end we first determine the relaxation eigenvalues
λ_k and eigenmodes _k,_k.
Note that we consider the eigensystem in the presence of a
reflecting wall. Without loss of generality we here explicitly treat
only the “left”
problem -∞<x≤ a, since the opposite “right” problem (denoted
later on with †) follows
by analogy.  As before we have the normalization ∫_-∞^a_k(x)_l(x) x=δ_kl.
Using {_k(x),_l(x)} we now determine the first passage eigenvalues μ_k as explained in the previous subsection.
The resulting first passage eigenvalues μ_k will
automatically contain only
first passage modes corresponding to the “left” problem. This procedure remarkably simplifies the numerical
determination of the first passage distribution, especially of those
modes that are faster than the slowest mode of the “right” problem
(i.e. absorption from the right), μ_k>μ_1^†, since a small number of modes M
in Eq. (<ref>) might otherwise be confused with `fantom' modes
from the opposite case †. In a related work <cit.> we
investigated the “left” first passage problem (case (I)) for a triple well potential in the presence of a reflecting
boundary, and found an excellent agreement between the analytical first
passage time distribution and computer simulations extending 
over many orders of magnitude in time using merely
M=40 relaxation modes.
Finally, we have to point out that solving an eigenvalue problem {λ_k,_k} with reflecting boundary condition is numerically easier than without reflecting boundary, i.e., natural boundaries albeit theoretically easier are numerically harder.


In the following section we apply these theoretical results to a discrete-state protein folding model and for the Ornstein-Uhlenbeck process.





§ EXAMPLES







 §.§ Discrete protein folding model



 






We consider a simple continuous-time Markov state model 
for a protein with three structural elements as shown in
Fig. <ref>. The protein starts from an initially unfolded state x_0=(0,0,0)≡I,
from which it is searching for the native state a=(1,1,1) through
intermediate states II-VII (see e.g., <cit.>).
Each arrow in Fig. <ref> indicates a possible transition x→ x' (x,x'=I,…,VIII, x≠ x')
that occurs with a Arrhenius type rate L(x',x)≡exp(F_x-B_xx'), where
F_x denotes the free energy of state x and B_xx'=B_x'x the
energy barrier along the transition link x↔ x'.
The resulting transition matrix satisfies detailed balance ln [L(x',x)/L(x,x')]=F_x-F_x' for all values of {F_x,B_xx'},
and naturally has negative diagonal elements
L(x,x)=-∑_x'≠ xL(x',x).

To test the power of the method from Sec. <ref>
we set up the 8×8 transition matrix
 with elements ⟨x'||x⟩=L(x',x) for a given set of energy barriers B_xx' and free energies F_x.
Then we carry out the eigendecomposition of , for which we first determine
the eigenvalues 0<λ_1,…, λ_7 (with λ_0=0)
corresponding to the zeros of the characteristic function
(<ref>), χ(-λ_k)=0. We then determine the right eigenvectors |_k⟩
by solving |_k⟩=-λ_k|_k⟩ for k=0,…,7.
The corresponding left eigenvectors, which solve ⟨_k|=-⟨_k|λ_k,
have components
⟨_k|x=𝒩^-1^F_x⟨x|_k,
where 𝒩=∑_x=I^VIII^F_x|⟨_k|x|^2 is a normalization factor. We take the function
F_k^* as defined in Eq. (<ref>), where
Ψ_k=⟨a|_k⟨_k|a and
μ̅_k=(λ_k+λ_k-1)/2 with k=0,…,7, and choose k^*(k)=k,k-1 according to Eq. (<ref>), which guarantees f(s,k)=F_k^*(k)(s)
to be negative at s=-μ̅_k.
The truncated Newton series (<ref>) involving  the first N terms
is then given by

    μ_k^N=μ̅_k+∑_n=1^N f_0(k)^n/f_1(k)^2n-1_n(k)/(n-1)!,

where _n(k)
is the determinant of the almost triangular matrix from Eq. (<ref>)
and 
f_i(k) is the ith derivative of f(s,k) at s=-μ,
with explicit formulas given in Eq. (<ref>).
The weights w_k(x_0) are determined using Eq. (<ref>),
i.e. by inserting μ_k→μ_k^N.
The calculations are performed for 100 randomly generated folding landscapes chosen 
as described in the caption to Fig. <ref>.







In Fig. <ref> we present the results for one particular realization of the folding landscape.
Fig. fig:folding_result1a displays the first passage time distribution
for N=6 and N=12 on a doubly-logarithmic scale. The solid line
represents the first passage time distribution obtained via a
numerical diagonalization of _a. The corresponding duality solutions nicely overlap with the numerical result even on relatively short time scales (see inset for a plot with linear scales).

Having obtained the full distribution of first passage times is
important for understanding kinetics in the so-called few encounter
limit <cit.>, in which for example 100 molecules are
simultaneously searching for a state a. This scenario is indeed
biologically relevant, for example, in the misfolding-triggered protein
aggregation, which in turn leads to numerous diseases (see
<cit.> for a more detailed discussion). Namely, as soon as the first protein
molecule spontaneously misfolds it creates a nucleation site for further
downhill misfolding and aggregation events ultimately leading to a
macroscopic insoluble toxic aggregate. 

In such a scenario the 
typical timescale of first arrivals will naturally be shifted towards
shorter timescales, thus requiring an accurate determination of the full
first passage statistics. 
Standard approaches focusing on the mean first passage time alone,
would therefore fail in the few encounter limit, whereas 
our new framework provides an accurate and consistent result (see also <cit.> and Fig. 2 therein for more details).




In Fig. <ref> we systematically analyze the deviation of the truncated
Newton series (<ref>) 
with respect to corresponding numerically obtained first passage eigenvalues μ_k
for 100 randomly generated
folding landscapes. For a given landscape the relative error is quantified in
terms of the dimensionless quantity ϵ≡|μ_k^N-μ_k|/μ_k, and Fig. fig:folding_result2a depicts the typical error
characterized by the median of the individual errors for all seven modes, respectively.
Note that for N=20 the relative error of the finite Newton series is typically below 10^-6.





Fig. fig:folding_result2b displays the maximal error out
of 100 randomly picked landscapes. We observe that larger errors can
occur if a first passage eigenvalue is located immediately after a gap in the relaxation spectrum.
The smaller error of the slowest first passage mode  μ_1
is due to the fact that μ_1 cannot be located after such a gap due to
the interlacing theorem (<ref>), which implies μ_1≤λ_1.
In this specific example the maximum relative error out of 
100 models randomly generated models
is found for the fifth mode (μ_5); the inset of Fig. fig:folding_result2b shows the corresponding P̃(a,s|a) as well as the result form the finite Newton
series with N=20 (see diamonds in the inset of Fig. fig:folding_result2b).
In this extreme scenario the weight of the fifth relaxation mode Ψ_5≪Ψ_l≠ 5
is almost negligible compared to other weights, leading to an almost
vanishing weight w_5, which would in turn require an increased number
of terms N entering the Newton series. Increasing the number of terms in the truncated Newton series from N=20 to N=200
reduces the deviation from  ϵ≃10^-1
to ϵ≃10^-2, the result for N=200 is marked by
the star in the inset of Fig. fig:folding_result2b.
Fig. <ref> readily demonstrates that the our
duality can be robustly and reliably applied to
all Markov state models.




 §.§ Ornstein-Uhlenbeck process


Let us now consider a linear Ornstein-Uhlenbeck, which corresponds to
a diffusion process in a harmonic potential β U(x)=ω
x^2/2. The corresponding Fokker-Planck operator reads =D_xω x+D_x^2.
The eigendecomposition of the relaxation process in the absence of
reflecting boundaries is well known. The respective eigenvalues
are given by λ_k=Dω k with the corresponding eigenfunctions  <cit.>

    _k(x)   ≡⟨x|_k=
      ^-ω x^2/2/√(2π/ω)H_k(x√(ω/2))/k!2^k,
    _k(x)   ≡⟨_k|x =H_k(x√(ω/2)),

where H_k is the kth Hermite polynomial. Although this process is extremely well studied, a closed-form analytical result for the first passage time distribution ℘_a(t|x_0) for any non-centered target position a≠ 0
remained elusive <cit.>.
We note that the well known analytical solution of the Laplace transform of the probability density ℘̃_a(s|x_0) in terms of Hermite polynomials <cit.> until now could only be analytically inverted to ℘_a(t|x_0) for the special case a=0 <cit.>.
Furtheremore, the exact large deviation limit ℘_a(t|x_0)≃ w_1(x_0)^-μ_1 t was just recently derived in <cit.>.
To obtain the full first passage time distribution we here use Eqs. (<ref>)-(<ref>) as follows.
Inserting Eq. (<ref>) into Eqs. (<ref>)-(<ref>)
yields the modifed diagonal of the propagator

    f(s,k)
        =(s+Dω k^*)∑_l=0^M
     ^-ω a^2/2/√(2π/ω)H_l(a√(ω/2))^2/l!2^l(s+Dω l)

where k^*(k)=k,k-1 is chosen according to Eq. (<ref>), which is equivalent to
f(-μ̅_k,k)<0 with  μ̅_k=Dω (k-1/2).
Note that we truncated the sum after M terms for the numerical evaluation, whereas
the exact formal result corresponds to M=∞.
The first line of Eq. (<ref>) is then simply given by f_0(k)=f(-μ̅_k,k) and the second line of Eq. (<ref>) becomes

    f_n≥ 1(k)=^-ω a^2/2/(Dω)^n√(2π/ω)∑_l=0
    l≠ k^*^M
      H_l(a√(ω/2))^2(l-k^*)/l!2^l(k-l-1/2)^n+1.

The kth first passage eigenvalue μ_k is determined by using the finite Newton series (<ref>), where the almost triangular matrix is taken from Eq. (<ref>), and the corresponding first passage weights w_k(x_0)
are determined using the residue theorem (<ref>). Note that our theory allows for the first time to determine analytically all first passage eigenvalues {μ_k} as well as the weights {w_k} and, therefore, also provides a complete solution to the first passage time density ℘_a(t|x_0).







Fig. <ref> depicts the results for the case, where the
absorbing point is set at a=2. Note that this scenario does not yet
correspond to the well-known high barrier Kramers regime. 
In Fig. fig:harmonica we compare the exact P̃(a,s|a) (solid blue line) with the finite approximation from Eq. (<ref>) using M=20 (dashed yellow line) and M=2000 (dash-dotted red line) relaxation modes, respectively.
The symbols represent the corresponding first passage modes (s=-μ_1,-μ_2,…).
Using only a small number of relaxation modes M=20 (see yellow
crosses) the zeros differ substantially from the respective
numerically obtained solution, which becomes, however, rather well approximated if we increase
the number of modes to M=2000 (see open red circles).
We note that such deviations of the first passage modes become particularly
inconvenient
for the modes that are marked by the vertical arrows in Fig. fig:harmonica. These first passage modes correspond to an absorption from the right, where
the corresponding weights  vanish w_k(x_0)=0 (here k=4,7,9,…) completely for all x_0≤ a, which, however, is only obtained in the limit M→∞.

This numerical truncation problem can be avoided completely if the relaxation process is considered with a 
reflecting boundary condition as explained in Sec. <ref>,
which automatically removes beforehand all zeros marked by the arrows in Fig. fig:harmonica (s=-μ_4,-μ_7,-μ_9,…).
Nevertheless, to illustrate the power and robustness of our duality approach we
proceed here without a reflecting wall and use M=2000.
In  Fig. fig:harmonicb we show the first passage time distribution
on a log-log scale (see inset for a linear scale) for three different starting positions x_0=0,0.5,1 (absorbing point a=2). The lines represent the first passage time distribution which is determined using our new method (with M=2000
relaxation modes) and the symbols represent the results ℘_a of a
numerical Laplace inversion of the renewal theorem (see figure caption for more details).
We find a perfect agreement between our new analytical method (lines)
and the numerical solution.
For comparison, we imposed a reflecting wall at the target in a related article <cit.> and obtained a similarly excellent agreement
between the duality solution and the simulated first passage time
density using a total of M=40 relaxation modes to quantify the first
passage time statistics for diffusion in a multi-well potential.
In either case, our new duality framework is exact for infinite M,
and hence the desired precision can be tuned at will.







§ CONCLUDING PERSPECTIVES


We rigorously established a duality between the relaxation and the
corresponding first passage processes
in terms of an interlacing of eigenvalues.  In other words, the time-scales at which a particle is
absorbed into the target are proven to interlace with the corresponding relaxation timescales. This duality allows us to understand first passage processes, both qualitatively and quantitatively, in terms of relaxation eigenmodes. For example,  spectral gaps in the relaxation spectrum
translate directly into spectral gaps in the first passage spectrum.
More explicitly, in effectively  one dimensional systems N gaps in the relaxation spectrum,
arising from N local (free) energy basins, translate into N-1 gaps in the first passage
time spectrum corresponding to the N-1 barriers separating the minima.
Most importantly, we established a duality that allows,  for the first time,
to determine exactly the first passage time distribution from
the corresponding relaxation spectrum.


Our theory is developed end tested
on both, continuous reversible Fokker-Planck dynamics and Markov state jump
processes in arbitrary dimensions.
For convenience and without loss of generality, we restricted the
applications of the duality for systems obeying Fokker-Planck dynamics
to effectively one dimensional problems.
An extension to more general models, for example, to diffusion on graphs would be straightforward, albeit rendering the calculations more cumbersome.


We tested and applied our theory to a discrete Markov state model 
of a simple protein folding landscape and the Ornstein-Uhlenbeck process, while
a continuous analogue of a  folding landscape are discussed elsewhere <cit.>.
Notably, we have derived, to the best of our knowledge, for the first time an exact and explicit analytical expression for the first passage time distribution of the Ornstein-Uhlenbeck process.

Looking forward it will be interesting and relevant to apply the
duality to the analysis of first passage processes on
graphs. Applications of the duality to narrow escape problems in
arbitrary dimensions <cit.>
will also be carried out in future studies. 

Finally, an extension of the framework to periodically or constantly driven
systems (i.e., irreversible Markovian dynamics), which goes beyond the long time limit that is presented
in <ref>,
will be particularly
challenging. Namely, there the interlacing theorem
cannot  be expected to hold anymore, since both eigenvalue spectra {λ_k} and {μ_k} can become complex valued.






§ EXPLICIT FORMULA FOR PRINCIPAL EIGENVALUE


In this appendix we simplify Eqs. (<ref>)-(<ref>) in the limit of a time-scale separation and for rare-event asymptotics for the principal
first passage eigenvalue μ_1.
We obtain a compact asymptotic expression of the principal first passage eigenvalue
μ̃_1≃μ_1, which is particular accurate if the time-scale of the slowest first passage eigenvalue is well separated from the time-scale of the slowest relaxation
mode (μ_1≪λ_1),
which inter alia
refines a previously proposed approximate link between the mean first passage time and the slowest relaxation mode <cit.>.

First, we redefine Eq. (<ref>) by setting 
k=1, k^*=0 and μ̅_1=0

    f(s)=P_eq(a)+∑_l≥ 1⟨a|_l⟨_l|as/s+λ_l,

where we dropped for convenience any argument with k since k=1
is assumed throughout this appendix.
The nth derivative of f at s=0
simplifies with Eq. (<ref>) to

    f_n=
     {[                         P_eq(a)                         if n=0,; n!(-1)^n+1∑_l≥1⟨a|_l⟨_l|a/λ_l^n                         if n≥1. ].

The almost triangular matrices Eq. (<ref>) become

    𝒜_n^i,j=f_i-j+2Θ(i-j+1)/(i-j+2)![n(i-j+1)Θ(j-2)
     +iΘ(1-j)+j-1],

where we have replaced f_n(k) by f_n from Eq. (<ref>).
Consequently, the Newton series (<ref>)
also simplifies to

    μ_1=∑_n=1^∞f_0^n/f_1^2n-1𝒜_n/(n-1)!.







If we now set f_3=f_4=…=0  in the almost triangular matrices (<ref>), that is 𝒜̃_̃ñ≡𝒜_n|_f_3=f_4=…=0,
the resulting matrix 𝒜̃_̃ñ becomes triangular,
implying that its determinant is simply given by the
product of the diagonal elements

    𝒜̃_̃ñ   =∏_i=1^n-1𝒜_n^i,i
    =(f_2/2)^n-1(2n-2)!/n!,

where we have inserted Eq. (<ref>) and evaluated the product
in the last step.
Replacing 𝒜→𝒜̃_̃ñ
in the
Newton series (<ref>) finally
yields exactly

    μ̃_1≡∑_n=1^∞f_0^n/f_1^2n-1𝒜̃_̃ñ/(n-1)!=f_1-√(f_1^2-2f_0f_2)/f_2.

Eq. (<ref>) is nothing but the
root of the second order Taylor expansion of f(s) around s=0 (i.e., the parabolic equation).
This approximation is quite accurate whenever μ̃_1≪λ_1.

If the target is located at a high energy barrier,
such that slowest first passage eigenvalue is exponentially suppressed by the 
(free) energy at the target (i.e.,  μ_1∝^-U(a)), Eq. (<ref>)
will lead to a quite accurate approximation μ̃_1, which can be seen in Fig. <ref>.
More precisely, in Fig. fig:muaa we depict the relative error |μ̃_1-μ_1|/μ_1
as function of the target-site energy E =a^2/2
for the Ornstein-Uhlenbeck process from Sec. <ref>
with U(x)=x^2/2 (ω=D=β=1).
Conversely, Fig. fig:muab displays the results the first passage time distribution
to four different target positions
a=1.5,2,2.5,3 for a particle starting from x_0=1.
The symbols represent histograms for ℘_a(t|x_0) deduced from 10^6 Brownian dynamics trajectories,
and the lines correspond to the large deviation asymptotic
℘_a(t|x_0)≃ w_1(x_0)μ̃_1^-μ̃_1 t deduced from (<ref>).
We conclude that the limit  μ_1≪λ_1 lead to both, a quite accurate approximation μ̃_1≃μ_1
and to an effectively single exponential decay μ_1≪μ_2 of the first passage statistics, which extends previous results <cit.> (see also <cit.>).



Moreover, a spectral gap such as λ_1≪λ_2 
will also render f_n≥ 2 from Eq. (<ref>) to be negligibly small
if the target is not located at the global minimum of the potential. For example,
a multi-barrier crossing, as the one studied in Ref. <cit.> (see Fig. 5 therein),
the principal first passage eigenvalue from Eq. (<ref>)
deviates less than two percent from the exact value μ_1, i.e.,
|μ_1-μ̃_1|/μ_1<0.02.
Notably, the approximation
(<ref>)
refines previous conjectures that the mean first passage time to escape from
the deepest potential basin corresponds to the first nonzero relaxation mode <cit.>.


We note that Eq. (<ref>) can be reformulated to give

    μ̃_1=σ_1/2σ_2[√(1+4P_eq(a)σ_2/σ_1^2)-1],

where we inserted f_0=P_eq(a) and defined σ_n=∑_l≥1⟨a|_l⟨_l|a/λ_l^n. This relation is equivalent to Eq. (17)
from a related article <cit.>.
 
Finally, we emphasize that Eqs. (<ref>) and (<ref>) apply also to  homogeneous irreversible Markov processes, i.e.,  the relations from this appendix are not restricted
to hold just for reversible Markov chains.





We thank Matteo Polettini for useful comments on our manuscript.
The financial support from the German Research
Foundation (DFG) through the Emmy Noether Program
“GO 2762/1-1” (to AG) is gratefully acknowledged.



§ REFERENCES


10

kram40
H. Kramers, “Brownian motion in a field of force and the diffusion model of
  chemical reactions,”
  http://dx.doi.org/10.1016/S0031-8914(40)90098-2Physica
  7 (1940) 284–304.

szab80
A. Szabo, K. Schulten, and Z. Schulten, “First passage time approach to
  diffusion controlled reactions,”
  http://dx.doi.org/10.1063/1.439715J. Chem. Phys. 72
  (1980) 4350–4357.

ben93
E. Ben-Naim, S. Redner, and F. Leyvraz, “Decay kinetics of ballistic
  annihilation,” http://dx.doi.org/10.1103/PhysRevLett.70.1890Phys. Rev. Lett. 70 (1993) 1890–1893.

osha95
G. Oshanin, A. Stemmer, S. Luding, and A. Blumen, “Smoluchowski approach for
  three-body reactions in one dimension,”
  http://dx.doi.org/10.1103/PhysRevE.52.5800Phys. Rev. E
  52 (1995) 5800–5805.

meji11
C. Mejía-Monasterio, G. Oshanin, and G. Schehr, “First passages for a search
  by a swarm of independent random searchers,”
  http://dx.doi.org/10.1088/1742-5468/2011/06/P06022J. Stat.
  Mech. (2011) P06022.

guer16
T. Guérin, N. Levernier, O. Bénichou, and R. Voituriez, “Mean
  first-passage times of non-Markovian random walkers in confinement,”
  http://dx.doi.org/10.1038/nature18272Nature 534
  (2016) 356–359.

li17
Y. Li, D. Debnath, P. K. Ghosh, and F. Marchesoni, “Nonlocality of relaxation
  rates in disordered landscapes,”
  http://dx.doi.org/10.1063/1.4976844J. Chem. Phys. 146 (2017) 084104.

haen90
P. Hänggi, P. Talkner, and M. Borkovec, “Reaction-rate theory: fifty years
  after Kramers,” http://dx.doi.org/10.1103/RevModPhys.62.251Rev.
  Mod. Phys. 62 (1990) 251–341.

redn01
S. Redner, http://dx.doi.org/10.1017/CBO9780511606014A guide to
  first-passage processes.
Cambridge University press, Cambridge, 2001.

metz14a
R. Metzler, G. Oshanin, and S. Redner, eds.,
  http://dx.doi.org/10.1142/9789814590297_fmatterFirst-Passage
  Phenomena and Their Applications.
World Scientific Publishing, Singapore, 2014.

beni14
O. Bénichou and R. Voituriez, “From first-passage times of random walks in
  confinement to geometry-controlled kinetics,”
  http://dx.doi.org/https://doi.org/10.1016/j.physrep.2014.02.003Phys. Rep. 539 (2014) 225–284.

kope88
R. Kopelman, “Fractal reaction kinetics,”
  http://dx.doi.org/10.1126/science.241.4873.1620Science
  241 (1988) 1620–1626.

avra00
D. ben Avraham and S. Havlin,
  http://dx.doi.org/10.1017/CBO9780511605826Diffusion and
  Reactions in Fractals and Disordered Systems.
Cambridge University Press, Cambridge, 2000.

rupp15
J.-F. Rupprecht, O. Bénichou, D. S. Grebenkov, and R. Voituriez, “Exit time
  distribution in spherically symmetric two-dimensional domains,”
  http://dx.doi.org/10.1007/s10955-014-1116-6J. Stat. Phys.
  158 (2015) 192–230.

greb16
D. S. Grebenkov, “Universal formula for the mean first passage time in planar
  domains,” http://dx.doi.org/10.1103/PhysRevLett.117.260201Phys.
  Rev. Lett. 117 (2016) 260201.

bres13
P. C. Bressloff and J. M. Newby, “Stochastic models of intracellular
  transport,” http://dx.doi.org/10.1103/RevModPhys.85.135Rev.
  Mod. Phys. 85 (2013) 135–196.

gode15
A. Godec and R. Metzler, “Optimization and universality of Brownian search in
  a basic model of quenched heterogeneous media,”
  http://dx.doi.org/10.1103/PhysRevE.91.052134Phys. Rev. E
  91 (2015) 052134.

vacc15
G. Vaccario, C. Antoine, and J. Talbot, “First-passage times in
  d-dimensional heterogeneous media,”
  http://dx.doi.org/10.1103/PhysRevLett.115.240601Phys. Rev.
  Lett. 115 (2015) 240601.

gode16a
A. Godec and R. Metzler, “First passage time distribution in heterogeneity
  controlled kinetics: going beyond the mean first passage time,”
  http://dx.doi.org/10.1038/srep20349Sci. Rep. 6
  (2016) 20349.

vile09
R. D. Vilela and B. Lindner, “Comparative study of different
  integrate-and-fire neurons: Spontaneous activity, dynamical response, and
  stimulus-induced correlation,”
  http://dx.doi.org/10.1103/PhysRevE.80.031909Phys. Rev. E
  80 (2009) 031909.

brau15
W. Braun, P. C. Matthews, and R. Thul, “First-passage times in
  integrate-and-fire neurons with stochastic thresholds,”
  http://dx.doi.org/10.1103/PhysRevE.91.052701Phys. Rev. E
  91 (2015) 052701.

bark14
E. Barkai, E. Aghion, and D. A. Kessler, “From the area under the Bessel
  excursion to anomalous diffusion of cold atoms,”
  http://dx.doi.org/10.1103/PhysRevX.4.021036Phys. Rev. X
  4 (2014) 021036.

sing06
A. Singer, Z. Schuss, and D. Holcman, “Narrow escape, part ii: The circular
  disk,” http://dx.doi.org/10.1007/s10955-005-8027-5J. Stat.
  Phys. 122 (2006) 465–489.

schu07
Z. Schuss, A. Singer, and D. Holcman, “The narrow escape problem for diffusion
  in cellular microdomains,”
  http://dx.doi.org/10.1073/pnas.0706599104Proc. Natl. Acad. Sci.
  USA 104 (2007) 16098–16103.

rein09
J. Reingruber and D. Holcman, “Gated narrow escape time for molecular
  signaling,” http://dx.doi.org/10.1103/PhysRevLett.103.148102Phys. Rev. Lett. 103 (2009) 148102.

pill10
S. Pillay, M. J. Ward, A. Peirce, and T. Kolokolnikov, “An asymptotic analysis
  of the mean first passage time for narrow escape problems: Part I:
  Two-dimensional domains,” http://dx.doi.org/10.1137/090752511Multiscale Model. Simul. 8 (2010) 803–835.

isaa16
S. A. Isaacson, A. J. Mauro, and J. Newby, “Uniform asymptotic approximation
  of diffusion to a small target: Generalized reaction models,”
  http://dx.doi.org/10.1103/PhysRevE.94.042414Phys. Rev. E
  94 (2016) 042414.

greb17
D. S. Grebenkov and G. Oshanin, “Diffusive escape through a narrow opening:
  new insights into a classic problem,”
  http://dx.doi.org/10.1039/C6CP06102HPhys. Chem. Chem. Phys.
  19 (2017) 2723–2739.

paly14
V. V. Palyulin, A. V. Chechkin, and R. Metzler, “Lévy flights do not always
  optimize random blind search for sparse targets,”
  http://dx.doi.org/10.1073/pnas.1320424111Proc. Natl. Acad. Sci.
  USA 111 (2014) 2931–2936.

gode17
A. Godec and R. Metzler, “First passage time statistics for two-channel
  diffusion,” http://dx.doi.org/10.1088/1751-8121/aa5204J. Phys.
  A: Math. Theor. 50 (2017) 084001.

beni11
O. Bénichou, C. Loverdo, M. Moreau, and R. Voituriez, “Intermittent search
  strategies,” http://dx.doi.org/10.1103/RevModPhys.83.81Rev.
  Mod. Phys. 83 (2011) 81–129.

maju01
S. N. Majumdar and A. J. Bray, “Spatial persistence of fluctuating
  interfaces,” http://dx.doi.org/10.1103/PhysRevLett.86.3700Phys.
  Rev. Lett. 86 (2001) 3700–3703.

maju02
S. N. Majumdar and A. Comtet, “Exact asymptotic results for persistence in the
  Sinai model with arbitrary drift,”
  http://dx.doi.org/10.1103/PhysRevE.66.061105Phys. Rev. E
  66 (2002) 061105.

bray13
A. J. Bray, S. N. Majumdar, and G. Schehr, “Persistence and first-passage
  properties in nonequilibrium systems,”
  http://dx.doi.org/10.1080/00018732.2013.803819Adv. Phys.
  62 (2013) 225–361.

jarz11
C. Jarzynski, “Equalities and inequalities: Irreversibility and the second law
  of thermodynamics at the nanoscale,”
  http://dx.doi.org/10.1146/annurev-conmatphys-062910-140506Annu.
  Rev. Condens. Matter Phys. 2 (2011) 329–351.

seif12
U. Seifert, “Stochastic thermodynamics, fluctuation theorems and molecular
  machines,” http://dx.doi.org/10.1088/0034-4885/75/12/126001Rep.
  Prog. Phys. 75 (2012) 126001.

broe15
C. Van den Broeck and M. Esposito, “Ensemble and trajectory thermodynamics: A
  brief introduction,”
  http://dx.doi.org/http://dx.doi.org/10.1016/j.physa.2014.04.035Physica A 418 (2015) 6–16.

neri17
I. Neri, É. Roldán, and F. Jülicher, “Statistics of infima and stopping
  times of entropy production and applications to active molecular processes,”
  http://dx.doi.org/10.1103/PhysRevX.7.011019Phys. Rev. X
  7 (2017) 011019.

fuch16
J. Fuchs, S. Goldt, and U. Seifert, “Stochastic thermodynamics of resetting,”
  http://dx.doi.org/10.1209/0295-5075/113/60009EPL 113 (2016) 60009.

rold17
É. Roldán and S. Gupta, “Path-integral formalism for stochastic resetting:
  Exactly solved examples and shortcuts to confinement,”
  http://dx.doi.org/10.1103/PhysRevE.96.022130Phys. Rev. E
  96 (2017) 022130.

garr17
J. P. Garrahan, “Simple bounds on fluctuations and uncertainty relations for
  first-passage times of counting observables,”
  http://dx.doi.org/10.1103/PhysRevE.95.032134Phys. Rev. E
  95 (2017) 032134.

ging17
T. R. Gingrich and J. M. Horowitz, “Fundamental bounds on first passage time
  fluctuations for currents,”
  http://dx.doi.org/10.1103/PhysRevLett.119.170601Phys. Rev.
  Lett. 119 (2017) 170601.

bara15
A. C. Barato and U. Seifert, “Thermodynamic uncertainty relation for
  biomolecular processes,”
  http://dx.doi.org/10.1103/PhysRevLett.114.158101Phys. Rev.
  Lett. 114 (2015) 158101.

kole07
G. Kolesov, Z. Wunderlich, O. N. Laikova, M. S. Gelfand, and L. A. Mirny, “How
  gene order is influenced by the biophysics of transcription regulation,”
  http://dx.doi.org/10.1073/pnas.0700672104Proc. Natl. Acad. Sci.
  USA 104 (2007) 13948–13953.

fras07
P. Fraser and W. Bickmore, “Nuclear organization of the genome and the
  potential for gene regulation,”
  http://dx.doi.org/10.1038/nature05916Nature 447
  (2007) 413–417.

beni10
O. Bénichou, C. Chevalier, J. Klafter, B. Meyer, and R. Voituriez,
  “Geometry-controlled kinetics,”
  http://dx.doi.org/10.1038/nchem.622Nat. Chem. 2
  (2010) 472–477.

meye11
B. Meyer, C. Chevalier, R. Voituriez, and O. Bénichou, “Universality classes
  of first-passage-time distribution in confined media,”
  http://dx.doi.org/10.1103/PhysRevE.83.051116Phys. Rev. E
  83 (2011) 051116.

gode16
A. Godec and R. Metzler, “Universal proximity effect in target search kinetics
  in the few-encounter limit,”
  http://dx.doi.org/10.1103/PhysRevX.6.041037Phys. Rev. X
  6 (2016) 041037.

bial12
W. Bialek, Biophysics: Searching for Principles.
Princeton University Press, 2012.

pulk13
O. Pulkkinen and R. Metzler, “Distance matters: The impact of gene proximity
  in bacterial gene regulation,”
  http://dx.doi.org/10.1103/PhysRevLett.110.198101Phys. Rev.
  Lett. 110 (2013) 198101.

dobs03
C. M. Dobson, “Protein folding and misfolding,”
  http://dx.doi.org/10.1038/nature02261Nature 426
  (2003) 884–890.

chit06
F. Chiti and C. M. Dobson, “Protein misfolding, functional amyloid, and human
  disease,”
  http://dx.doi.org/10.1146/annurev.biochem.75.101304.123901Annu.
  Rev. Biochem. 75 (2006) 333–366.

zhen13
W. Zheng, N. P. Schafer, and P. G. Wolynes, “Free energy landscapes for
  initiation and branching of protein aggregation,”
  http://dx.doi.org/10.1073/pnas.1320483110Proc. Natl. Acad. Sci.
  USA 110 (2013) 20515–20520.

yu15
H. Yu, D. R. Dee, X. Liu, A. M. Brigley, I. Sosova, and M. T. Woodside,
  “Protein misfolding occurs by slow diffusion across multiple barriers in a
  rough energy landscape,”
  http://dx.doi.org/10.1073/pnas.1419197112Proc. Natl. Acad. Sci.
  USA 112 (2015) 8308–8313.

dee16
D. R. Dee and M. T. Woodside, “Comparing the energy landscapes for native
  folding and aggregation of PrP,”
  http://dx.doi.org/10.1080/19336896.2016.1173297Prion 10 (2016) 207–220.

nguy17
B. Nguyen, D. Hartich, U. Seifert, and P. De Los Rios, “Thermodynamic bounds
  on the ultra- and infra-affinity of Hsp70 for its substrates,”
  http://dx.doi.org/10.1016/j.bpj.2017.06.010Biophys. J.
  113 (2017) 362–370.

hart18_arxiv
D. Hartich and A. Godec, “Duality between relaxation and first passage in
  reversible Markov dynamics: Rugged energy landscapes disentangled,” https://dx.doi.org/10.1088/1367-2630/aaf038New J. Phys. 20 (2018) 112002.

schn17a
D. Schnoerr, G. Sanguinetti, and R. Grima, “Approximation and inference
  methods for stochastic biochemical kinetics—a tutorial review,”
  http://dx.doi.org/10.1088/1751-8121/aa54d9J. Phys. A: Math.
  Theor. 50 (2017) 093001.

muns09
B. Munsky, I. Nemenman, and G. Bel, “Specificity and completion time
  distributions of biochemical processes,”
  http://dx.doi.org/10.1063/1.3274803J. Chem. Phys. 131 (2009) 235103.

bel10
G. Bel, B. Munsky, and I. Nemenman, “The simplicity of completion time
  distributions for common complex biochemical processes,”
  http://dx.doi.org/10.1088/1478-3975/7/1/016003Phys. Biol.
  7 (2010) 016003.

grim17
R. Grima and A. Leier, “Exact product formation rates for stochastic enzyme
  kinetics,” http://dx.doi.org/10.1021/acs.jpcb.6b08891J. Phys.
  Chem. B 121 (2017) 13–23.

schn17
D. Schnoerr, B. Cseke, R. Grima, and G. Sanguinetti, “Efficient low-order
  approximation of first-passage time distributions,”
  http://dx.doi.org/10.1103/PhysRevLett.119.210601Phys. Rev.
  Lett. 119 (2017) 210601.

webe17
M. F. Weber and E. Frey, “Master equations and the theory of stochastic path
  integrals,” Rep. Prog. Phys. 80 (2017) 046601.

biro01
G. Biroli and J. Kurchan, “Metastable states in glassy systems,”
  http://dx.doi.org/10.1103/PhysRevE.64.016101Phys. Rev. E
  64 (2001) 016101.

tana03
S. Tănase-Nicola and J. Kurchan, “Topological methods for searching
  barriers and reaction paths,”
  http://dx.doi.org/10.1103/PhysRevLett.91.188302Phys. Rev. Lett.
  91 (2003) 188302.

tana04
S. Tănase-Nicola and J. Kurchan, “Metastable states, transitions, basins
  and borders at finite temperatures,”
  http://dx.doi.org/10.1023/B:JOSS.0000041739.53068.6aJ. Stat.
  Phys. 116 (2004) 1201–1245.

schu79
Z. Schuss and B. J. Matkowsky, “The exit problem: A new approach to diffusion
  across potential barriers,” http://dx.doi.org/10.1137/0136043SIAM J. Appl. Math. 36 (1979) 604–623.

matk81
B. J. Matkowsky and Z. Schuss, “Eigenvalues of the Fokker–Planck operator
  and the approach to equilibrium for diffusions in potential fields,”
  http://dx.doi.org/10.1137/0140020SIAM J. Appl. Math. 40 (1981) 242–254.

gron08
R. Grone, K. H. Hoffmann, and P. Salamon, “An interlacing theorem for
  reversible Markov chains,”
  http://dx.doi.org/10.1088/1751-8113/41/21/212002J. Phys. A:
  Math. Theor. 41 (2008) 212002.  
  
  
kamp07
N. G. van Kampen,
  http://dx.doi.org/10.1016/B978-044452965-7/50006-4Stochastic
  Processes in Physics and Chemistry.
North-Holland Personal Library. Elsevier, Amsterdam, 3 ed., 2007.

kada68
L. P. Kadanoff and J. Swift, “Transport coefficients near the critical point:
  A master-equation approach,”
  http://dx.doi.org/10.1103/PhysRev.165.310Phys. Rev. 165 (1968) 310–322.

gard04
C. W. Gardiner, http://dx.doi.org/10.1007/978-3-662-05389-8Handbook of Stochastic Methods.
Springer, Berlin, 3 ed., 2004.

sieg51
A. J. F. Siegert, “On the first passage time probability problem,”
  http://dx.doi.org/10.1103/PhysRev.81.617Phys. Rev. 81 (1951) 617–623.

pole17
M. Polettini and M. Esposito, “Effective thermodynamics for a marginal
  observer,” http://dx.doi.org/10.1103/PhysRevLett.119.240601Phys. Rev. Lett. 119 (2017) 240601.

bisk17
G. Bisker, M. Polettini, T. R. Gingrich, and J. M. Horowitz, “Hierarchical
  bounds on entropy production inferred from partial information,”
  http://dx.doi.org/10.1088/1742-5468/aa8c0dJ. Stat. Mech. (2017)
  093210.



keil64
J. Keilson, “A review of transient behavior in regular diffusion and
  birth-death processes,” http://dx.doi.org/10.2307/3211858J.
  Appl. Prob. 1 (1964) 247–266.

cahi02
N. D. Cahill, J. R. D'Errico, D. A. Narayan, and J. Y. Narayan, “Fibonacci
  determinants,” http://dx.doi.org/10.2307/1559033The College
  Mathematics Journal 33 (2002) 221–225.

jia18
J.-T. Jia, “Numerical algorithms for the determinant evaluation of general
  Hessenberg matrices,”
  http://dx.doi.org/10.1007/s10910-017-0794-0J. Math. Chem.
  56 (2018) 247–256.

titc62
E. C. Titchmarsh, Eigenfunction expansions associated with second-order
  differential equation.
Oxford University Press, London, 2 ed., 1962.

prin11
J.-H. Prinz, B. Keller, and F. Noé, “Probing molecular kinetics with Markov
  models: metastable states, transition pathways and spectroscopic
  observables,” http://dx.doi.org/10.1039/C1CP21258CPhys. Chem.
  Chem. Phys. 13 (2011) 16912–16927.

bowm14
G. R. Bowman, V. S. Pande, and F. Noé, eds.,
  http://dx.doi.org/10.1007/978-94-007-7606-7An Introduction to
  Markov State Models and Their Application to Long Timescale Molecular
  Simulation.
Springer Netherlands, Dordrecht, 2014.

alil05
L. Alili, P. Patie, and J. L. Pedersen, “Representations of the first hitting
  time density of an Ornstein-Uhlenbeck process,”
  http://dx.doi.org/10.1080/15326340500294702Stochastic Models
  21 (2005) 967–980.

greb15
D. S. Grebenkov, “First exit times of harmonically trapped particles: a
  didactic review,”
  http://dx.doi.org/10.1088/1751-8113/48/1/013001J. Phys. A: Math.
  Theor. 48 (2015) 013001.

nybe16
M. Nyberg, T. Ambjörnsson, and L. Lizana, “A simple method to calculate
  first-passage time densities with arbitrary initial conditions,”
  http://dx.doi.org/10.1088/1367-2630/18/6/063019New J. Phys.
  18 (2016) 063019.
  
  
darl53
D. A. Darling and A. J. F. Siegert, “The first passage problem for a continuous Markov process,”
  https://www.jstor.org/stable/2236785Ann. Math. Stat.
  24 (1953) 624–639.

verg18
M. Vergassola, V. E. Deneke, and S. Di Talia “Mitotic waves in the early embryogenesis of Drosophila: Bistability traded for speed,”
  http://dx.doi.org/10.1073/pnas.1714873115Proc. Natl. Acad. Sci. USA 115 (2018) E2165.
  
  
bere04
A. Berezhkovskii and A. Szabo, “Ensemble of transition states for two-state
  protein folding from the eigenvectors of rate matrices,”
  http://dx.doi.org/10.1063/1.1802674J. Chem. Phys. 121 (2004) 9186.

bico00
D. J. Bicout and A. Szabo, “Entropic barriers, transition states, funnels, and
  exponential protein folding kinetics: A simple model,”
  http://dx.doi.org/10.1110/ps.9.3.452Protein Sci. 9
  (2000) 452–465.





 