Stochastic path power and the Laplace transform

Functional- or path-integral techniques were first developed by Wiener in the 1920s as a means to study Brownian motion. Since Feynman's seminal work in the 1950s however [1], they have been more widely known for their use in quantum mechanics and quantum field theory. They have also been used to great effect in classical statistical physics and field theory; see for example [2, 3]. Their application to the modern study of classical stochastic processes can be traced back to the work of Onsager and Machlup [4], Stratonovich [5], and Graham [6]. Their asymptotic evaluation in the weak noise limit provides a powerful method for the study of rare-event processes such as barrier crossing rates. In this approximation, the path integral is analogous to the ray (or eikonal, or WKB) method for the corresponding Smoluchowski equation [7]. The generalization to correlated noise was accomplished in the 1990s by McKane and co-workers [8, 9], and a general field-theory based formalism was developed by Hochberg et al [10]. For an introduction to path integral methods for stochastic processes, see Wio [11], and for path integral methods in general Schulman [12] and Kleinert [13]. A more mathematical treatment in the context of large deviations theory was pioneered by Friedlin and Wentzell [14].

For the quantum-mechanical case, in the semiclassical limit of small $\hbar$, the dominant system trajectories are small fluctuations about the classical paths. In this limit, even the wholly quantum-mechanical phenomenon of barrier penetration can be understood through classical trajectories in the inverted potential $V\to-V$, known as instantons [15]. An analogous correspondence also exists for stochastic processes, where most-probable paths in a potential V correspond to Hamiltonian trajectories in an effective potential $-|\nabla V|^2$ [16, 17] and much of the formalism can be carried over [18]. Action-minimization methods to find the most probable path through higher-dimensional potentials have been developed [19–24]. The majority of previous work has focused on determining the infinite-time paths, and hence average transition rates. However, the full temporal information need not necessarily be discarded, and indeed is required for the calculation of quantities such as first-passage densities, and for situations where local equilibrium cannot be assumed.

In this paper, after a brief review of the formalism, we re-examine the stochastic–deterministic correspondence, and identify the conserved quantity in the effective mechanics as the Laplace parameter in the time-transformed diffusive dynamics. We illustrate the utility of the approach with the simple examples of a particle confined to an interval in one dimension, and the harmonic oscillator.

The starting point is an overdamped Markovian stochastic process driven by Gaussian white noise, defined by the Langevin equation (see e.g. [25])

$$\begin{align} \gamma\dot x\left(t\right) = -V^{^{\prime}}\left(x\right) + \xi\left(t\right);\;\;\langle\xi\left(t\right)\xi\left(t^{^{\prime}}\right)\rangle = 2D\delta\left(t-t^{^{\prime}}\right) \end{align} \tag{ 1 }$$

or

$$\begin{align} \gamma dX_t = -V^{^{\prime}}\left(X_t\right)dt+\sqrt{2D}dW, \end{align} \tag{ 2 }$$

with dW the increments of the standard Wiener process. V(x) is the potential, D is the noise strength, and γ is a dissipation parameter. Here γ represents the friction rather than the mobility, so it multiplies the $\dot x$ term rather than the force term. The fluctuation-dissipation theorem gives $D = \gamma k_\textrm{B}T$. The first equation above is not valid in the limit of zero correlation time, but it suggests how to define an action based on the noise probability density functional ${\mathcal P}[\xi]$:

$$\begin{align} {\mathcal P}\left[\xi\right]\sim\exp\left( -\frac1{4D}\int_0^t \xi\left(t^{^{\prime}}\right)^2{\mathrm{d}} t^{^{\prime}}\right). \end{align} \tag{ 3 }$$

This in turn leads to the following expression for the weight ${\mathcal P}[x]$ attached to a stochastic path x(t):

$$\begin{align} {\mathcal P}\left[x\right]\sim\exp\left( -\frac1{4D}\int_0^t \left(\gamma\dot x+V^{^{\prime}} \right)^2{\mathrm{d}} t^{^{\prime}}\right)\equiv \exp\left( -\frac{{\mathcal S}\left[x\right]}{4D}\right), \end{align} \tag{ 4 }$$

defining the stochastic action ${\mathcal S}$. This is only valid for delta-correlated white noise. Gaussian coloured noise, with a finite correlation time, has two time integrals in the exponent of equation (3), and hence a non-local stochastic action [8, 10]. In this work we keep D constant, but a coordinate-dependent D could be introduced by changing variables to a system with constant noise strength and a modified potential [26].

The probability density function $\rho(x,t)$ satisfies the Fokker–Planck–Smoluchowski equation

$$\begin{align} \frac{\partial \rho}{\partial t} = \frac{\partial}{\partial x}\left(\frac{\rho V^{^{\prime}}}{\gamma}+D\frac{\partial \rho}{\partial x}\right);\;\;\rho\left(x,0\right) = \delta\left(x-x_0\right), \end{align} \tag{ 5 }$$

and $\rho(x,t)$ corresponds to the conditional transition probability density $P(x,t|x_0,0)$. When the initial density is a delta function, this is the Green's function for the Smoluchowski equation. Many approaches to solving this equation have been developed; see [26] for a thorough treatment. The approach we follow here has most in common with the WKB method (appendix B), as it exploits the smallness of the noise strength D. The density can be written in terms of a functional integral:

$$\begin{align} \rho\left(x,t\right) = P\left(x,t|x_0,0\right) = \int{\mathcal D} x\,{\mathcal J}\left[x\right] \exp\left( -\frac{{\mathcal S}\left[x\right]}{4D}\right), \end{align} \tag{ 6 }$$

where ${\mathcal J} = \exp\left( \frac1{2}\gamma\int_0^tV^{^{\prime\prime}}(x(t^{^{\prime}}))\textrm{d} t^{^{\prime}}\right)$ is the functional Jacobian arising from the implicit change of variables $\xi\to x$. This follows from the Stratonovich interpretation of (1), see e.g. [11, 25]. If we interpreted (1) as an Itô SDE, this Jacobian would equal one. However, the cross term in the action integral would then be an Itô integral:

$$\begin{align} \int_0^t 2\gamma V^{^{\prime}}\dot x\,{\mathrm{d}} t^{^{\prime}} = 2\gamma \int_{x_0}^x V^{^{\prime}}\,{\mathrm{d}} x^{^{\prime}} = 2\gamma\left(V\left(x\right)-V\left(x_0\right)\right)+2\gamma D\int_0^tV^{^{\prime\prime}}\left(x\left(t^{^{\prime}}\right)\right)\,{\mathrm{d}} t^{^{\prime}}, \end{align} \tag{ 7 }$$

resulting in an additional term exactly equivalent to the non-unit Jacobian, as given above. A detailed discussion of this in the context of constructing a path-integral calculus is given in [27]. The integral is taken over paths $x(t^{^{\prime}})$ satisfying $x(0) = x_0,x(t) = x$, and an infinite normalization constant (arising from the '∼' in (3) above) has been absorbed into the measure. As D → 0, the integral is dominated by paths lying close to the smooth minimizer of ${\mathcal S}$, which we call $x_*(t)$. This satisfies the Euler–Lagrange equation for ${\mathcal S}$:

$$\begin{align} 2\gamma^2\ddot x_* = 2V^{^{\prime}}V^{^{\prime\prime}}, \end{align} \tag{ 8 }$$

that is, the conservative Hamiltonian dynamics of a particle of mass $2\gamma^2$ moving in the effective potential

$$\begin{align} F = -V^{^{\prime}}{}^2. \end{align} \tag{ 9 }$$

This remarkable correspondence, whereby the most probable trajectories of an overdamped stochastic process can be mapped onto an effective energy-conserving Hamiltonian system, was noted in [18] and first discussed in detail in [17].

This dynamics has a conserved quantity

$$\begin{align} H = \gamma\dot x_*^2 - V^{^{\prime}}\left(x_*\right)^2/\gamma\longrightarrow \dot x_*^2 - V^{^{\prime}} {}^2, \end{align} \tag{ 10 }$$

($\gamma^{-1}\times$) the energy of the effective system, which has dimensions of power. We set γ = 1 and call H the path power. If we had used the Itô interpretation with ${\mathcal J} = 1$, the effective potential would be $F = -V^{^{\prime}} {}^2 +2DV^{^{\prime\prime}}$. Both choices have been employed in the literature cited above, and whilst they are formally equivalent (as they must be for a constant-strength additive noise process), the dominant paths as D → 0 that we discuss below would be slightly different in cases where $V^{^{\prime\prime}}\neq 0$. Using $F = -V^{^{\prime}} {}^2$ facilitates a more straightforward interpretation of the conserved quantity H, which we will discuss in the next section. Also, it keeps the 'classical' action independent of the noise strength, and so separates the orders in D in analogy with the quantum-mechanical version.

As D → 0, we can approximate the expression for ρ by restricting the integral in (6) to paths 'near' the minimizer $x_*$. More precisely, quadratic fluctuations around $x_*$ may be integrated over by writing a general path x(t) as $x(t) = x_*(t)+y(t)$ and expanding the action to second order in y(t). This results in a Gaussian functional integral that can be performed explicitly, and yields the following expression for ρ:

$$\begin{align} \rho\left(x,t\right) = \left(4\pi D\,{\mathrm{det}}\,\hat L \right)^{-1/2}{\mathcal J}\left[x_*\right]\exp\left( -\frac{{\mathcal S}\left[x_*\right]}{4D}\right), \end{align} \tag{ 11 }$$

valid as D → 0. The operator $\hat L$ is the second variation of the action, and its determinant (the 'fluctuation determinant') can be evaluated by various techniques; see for example [12, 28] and appendix A. If more than one solution for $x_*$ exists, ρ will contain a term similar to the above for each of them.

According to equation (10), when $\dot x_* = - V^{^{\prime}}$, H = 0. This is the deterministic path, i.e. the solution to the noiseless equation of motion. It corresponds to infinite-time relaxation towards equilibrium, requires no fluctuations, and has action zero. H is also zero for $\dot x_* = + V^{^{\prime}}$, which is the most probable hill-climbing path in the long-time limit (see next section). Minimizing the action conditioned on a finite time results in a finite-H smooth minimizer $x_*$. In [17], paths with H ≠ 0 were interpreted as having a constant excess kinetic energy, added or subtracted along the path. Such smooth paths have measure zero in the space of non-differentiable stochastic paths over which we are integrating (almost all of which have H ≠ 0), but as D → 0, the stochastic paths satisfying the boundary conditions $x(0) = x_0,x(t) = x$ concentrate in a tube around $x_*$ [5]. Integrating over the quadratic fluctuations around $x_*$ corresponds to this tube, and approximates the path integral as D → 0.

Inserting the minimizer $x_*$ into the action functional ${\mathcal S}[x]$, and using $H = \dot x_*^2-V^{^{\prime}}(x_*)^2$ gives Hamilton's principal function $S(x,t)$ [17]:

$$\begin{eqnarray} S\left(x,t\right) = {\mathcal S}\left[x_*\left(t\right)\right] & = & \int_0^t \left(\dot x_* + V^{^{\prime}}\right)^2{\mathrm{d}} t^{^{\prime}} = 2\Delta\! V +\int_0^t \left(\dot x_*^2 + V^{^{\prime}} {}^2\right){\mathrm{d}} t^{^{\prime}} \nonumber\\ & = & 2\Delta\! V + \int_0^t \left(2\dot x_*^2 - H\right){\mathrm{d}} t^{^{\prime}} = 2\Delta\! V - Ht + 2\int_0^t \dot x_* \frac{{\mathrm{d}} x_*}{{\mathrm{d}} t^{^{\prime}}}{\mathrm{d}} t^{^{\prime}} \nonumber\\ & = & 2\Delta\! V - Ht + 2\int_{p}\sqrt{H + V^{^{\prime}} {}^2}\,{\mathrm{d}} s. \end{eqnarray} \tag{ 12 }$$

The term $\Delta\! V = V(x) - V(x_0)$ is path-independent and results from the total derivative in the action. s is the arc length along the path p, satisfying equation (8) — for a simple direct path in one dimension, this is simply from x₀ to x₁ along the x-axis [17]; for paths involving turns, or in higher dimensions, it is more complicated. H is as above, and is defined implicitly via

$$\begin{align} 0 = \frac{\partial S}{\partial H} = -t + \int_{p}\frac{{\mathrm{d}} s}{\sqrt{H + V^{^{\prime}} {}^2}}, \end{align} \tag{ 13 }$$

as can also be seen by integrating the energy equation of the effective classical mechanics $\dot x_* = \pm\sqrt{H+V^{^{\prime}}(x_*)^2}$. The fluctuation determinant det$\,\hat L$ is given by (more details are given in appendix A)

$$\begin{eqnarray} {\mathrm{det}}\,\hat L & = & \sqrt{H+V^{^{\prime}}\left(x_0\right)^2}\sqrt{H+V^{^{\prime}}\left(x\right)^2}\; \int_{p}\frac{{\mathrm{d}} s}{\left(H + V^{^{\prime}} {}^2\right)^{3/2}}\nonumber\\ & = & \sqrt{H+V^{^{\prime}}\left(x_0\right)^2}\sqrt{H+V^{^{\prime}}\left(x\right)^2} \; \, 2\,\left|\frac{\partial^2 S}{\partial H^2}\right|. \end{eqnarray} \tag{ 14 }$$

The appearance of the derivatives of S with respect to H is highly suggestive, and leads to our main result. Writing S as a Legendre transform

$$\begin{eqnarray} S\left(x,t\right) & = & -Ht + W\left(x,H\right), \nonumber\\ W\left(x,H\right) & = & 2\Delta\! V+ 2\int_{p}\sqrt{H + V^{^{\prime}} {}^2}\,{\mathrm{d}} s \end{eqnarray} \tag{ 15 }$$

defines W (Hamilton's characteristic function for the effective mechanics). This can be identified with the Hamilton–Jacobi equation of the effective classical mechanics, which was derived in [17], but when exponentiated resembles the integrand of an inverse Laplace transform from H to t. The Legendre transform exchanges the dependence on the conjugate variables t and H, which for the extremal path $x_*$ are related by equation (13). At fixed t, $\partial^2S/\partial H^2 = \partial^2W/\partial H^2$, and we can rewrite the D → 0 approximation for the density in equation (11) as

$$\begin{align} \rho\left(x,t\right) = A \exp\left( +\frac{Ht}{4D} \right)\exp\left( -\frac{W\left(x,H\right)}{4D} \right)\sqrt{\frac{8\pi D}{|W^{^{\prime\prime}}|}}, \end{align} \tag{ 16 }$$

where H satisfies (10) and the prefactor $A(x,H) = {\mathcal J}/\sqrt{\sqrt{H+V^{^{\prime}}(x_0)^2}\sqrt{H+V^{^{\prime}}(x)^2}}$. Now consider

$$\begin{align} \bar \rho\left(x;\frac{h}{4D}\right) = A\left(x,h\right) \exp\left( - \frac{W\left(x,h\right)}{4D} \right), \end{align} \tag{ 17 }$$

with $h\in \mathbb{C}$, as a candidate for the Laplace-transformed density. The time-domain density would then be given by

$$\begin{align} \rho\left(x,t\right) = \frac1{2\pi i}\int_{\mathrm{Br}} \exp\left( +\frac{ht}{4D} \right)A\left(x,h\right) \exp\left( -\frac{W\left(x,h\right)}{4D} \right) \frac{{\mathrm{d}} h}{4D}. \end{align} \tag{ 18 }$$

Taking the Bromwich contour as $h = H + iy;H\in\mathbb{R}^+,y\in(-\infty,\infty)$ and performing the integral by steepest descents as D → 0 gives equation (16). So, we can identify $H/4D$ as the real part of the Laplace parameter within the steepest descents approximation. This can be also be recovered from a WKB analysis of the Laplace transform of (5) (see appendix B).

The identification of the (noise-scaled) path power as the Laplace parameter, the variable conjugate to time, leads to a number of observations. Firstly, by the final value theorem, the H → 0 limit corresponds to $t\to\infty$. For an uphill path segment, where $\sqrt{V^{^{\prime}} {}^2} = +V^{^{\prime}}$, $W\to 4\Delta\! V$, and the familiar Kramers form $\exp\left( -\Delta\! V/D\right)$ emerges—the long time average rate for a process driven by noise of strength D to climb a potential barrier of height $\Delta\! V$. This relies on the assumption of quasi-equilibrium being reached at the bottom of the barrier, which we identify with H → 0. However, all the finite-time, non-equilibrium information remains when H is left finite. The awkwardness of the implicit definition of H(t) in (10) is no longer an issue when working in the Laplace domain. This is particularly helpful when multiple trajectories' contributions need to be summed, as the following example demonstrates.

Consider a particle diffusing on the x-axis, with $x = x_0\in (a,b)$ at t = 0, and reflecting boundaries at $x = a,b$ with b > a. This situation is analogous to the 'particle in a box' problem from elementary quantum mechanics. The density $\rho(x,t)$ satisfies the diffusion equation with zero-flux boundary conditions:

$$\begin{align} \frac{\partial \rho}{\partial t} = D\frac{\partial ^2\rho}{\partial x^2};\;\;\rho\left(x,0\right) = \delta\left(x-x_0\right);\;\;\frac{\partial \rho}{\partial x}\left(a,t\right) = 0 = \frac{\partial \rho}{\partial x}\left(b,t\right), \;\text{or} \end{align} \tag{ 19 }$$

$$\begin{equation} s\bar\rho - \delta\left(x-x_0\right) = D\frac{\partial ^2\bar\rho}{\partial x^2} \end{equation} \tag{ 20 }$$

for the transformed density $\bar\rho(x;s)$, with Laplace parameter s identified as $H/4D$. Since V = 0 here, the effective potential is also 0, the Jacobian ${\mathcal J} = 1$, and everything can be calculated explicitly. $A = H^{-1/2}$, and $W = 2\sqrt H |p|$ where $|p|$ is the total distance travelled, i.e. $|x-x_0|$ for a simple direct path from x₀ to x. Such a path gives

$$\begin{align} \bar\rho\left(x;s\right) = \frac{{\mathrm{e}}^{-2\sqrt{H}|x-x_0|/4D}}{\sqrt{H}} = \frac{{\mathrm{e}}^{-|x-x_0|\sqrt{s/D}}}{\sqrt{4Ds}} = {\mathrm{L.T.}}\; \left( \frac{{\mathrm{e}}^{-\left(x-x_0\right)^2/4Dt}}{\sqrt{4\pi Dt}}\right) \end{align} \tag{ 21 }$$

as expected (this can also be obtained directly from the untransformed expressions in the previous section). So, a simple direct path from x₀ to x gives the expected solution for the infinite region $a\to -\infty,b\to\infty$. With reflecting boundaries at $x = a,b$, however, there are four Hamiltonian trajectories that reach reach x from x₀: a bounce off x = b (which we call 'type 1'), a bounce off x = a ('type 2'), and a bounce off both ('type 3'), plus the simple direct path ('type 0'), as shown in figure 1. These have actions, numbered by type, given by

$$\begin{eqnarray} W_0/2\sqrt H & = & |x-x_0|\to x_0-a,\; b-x_0\;\;{\mathrm{as}}\;\;x\to a,b\nonumber\\ W_1/2\sqrt H & = & 2b-x_0-x\to 2b-a-x_0,b-x_0\;\;{\mathrm{as}}\;\;x\to a,b\nonumber\\ W_2/2\sqrt H & = & x+x_0-2a\to x_0-a,b-2a+x_0\;\;{\mathrm{as}}\;\;x\to a,b\nonumber\\ W_3/2\sqrt H & = & |x_0-x|+2b-2a\to 2b-a-x_0,x_0+b-2a\;\;{\mathrm{as}}\;\;x\to a,b \end{eqnarray} \tag{ 22 }$$

where we also give the limiting values of the actions at the boundaries. This leads to

$$\begin{align} \bar\rho = H^{-1/2}\left(\exp \left(-W_0/4D\right) + \exp \left(-W_1/4D\right) + \exp \left(-W_2/4D \right) + \exp\left( -W_3/4D\right) \right). \end{align} \tag{ 23 }$$

For x > 0, types 0 and 2 are increasing functions of x while types 1 and 3 are decreasing. This simply corresponds to the two independent solutions of the second-order ODE, and introduces relative sign changes in the x-derivatives. By linearity, $\bar\rho$ as above satisfies the transformed diffusion equation, and is zero-flux at x = a or b because of the way the terms in $\partial\bar\rho/\partial x$ cancel in pairs: $W_1(a) = W_2(a); W_3(a) = W_4(a)$ and $W_1(b) = W_3(b);W_2(b) = W_4(b)$. The initial condition can be checked (without Laplace inversion) using the initial value theorem:

$$\begin{align} \lim_{s\to\infty} s\bar\rho\left(x,s\right) = \lim_{t\to 0}\rho\left(x,t\right). \end{align} \tag{ 24 }$$

Because $W_0\lt W_{1,2,3}$ at fixed H, as s (or H) goes to infinity, only the first term contributes to the limit, recovering the infinite region solution, and hence the short-time limit is correctly captured. However, the long-time limit is incorrect: the final value theorem

$$\begin{align} \lim_{s\to 0} s\bar\rho\left(x,s\right) = \lim_{t\to \infty}\rho\left(x,t\right) \end{align} \tag{ 25 }$$

gives

$$\begin{align} \frac{H}{4D}\,\bar\rho\sim\frac{\sqrt{H}}{D}\; \text{as}\; H\to 0;\; \rho\left(x,t\to\infty\right)\to 0, \end{align} \tag{ 26 }$$

as in the infinite region case, whereas the reflecting boundaries should give a uniform distribution in this limit: $\rho(x,t\to\infty) = 1/(b-a)$.

Zoom In Zoom Out Reset image size

This apparent contradiction can be resolved by including the contributions of more paths. For each of the types of path above, we can add n whole circuits from x₀ to x₀. Each of these will accrue additional action $n\times 2\sqrt H\times 2(b-a)$. In the transformed domain, we work with a fixed H, so all these paths have the same H (in contrast to them all having the same t in the time domain). This means that the exponentials form a geometric series, and the multi-bounce paths can be summed over explicitly, resulting in

$$\begin{eqnarray} & \bar \rho = \frac{\exp\left( -W_0/4D\right) + \exp\left( -W_1/4D\right) + \exp\left( -W_2/4D \right) + \exp\left( -W_3/4D \right)}{\sqrt H\left(1 - \exp -\sqrt H \left(b-a\right)/D\right)},\nonumber\\ & \text{with}\;\; \lim_{H\to 0}\, \frac{H}{4D}\,\bar\rho \; = \frac1{b-a} = \rho\left(x,t\to\infty\right)\;\;\text{as required.} \end{eqnarray} \tag{ 27 }$$

Whilst the weak noise approximation inherent in the saddle-point evaluation of path integrals is exact in the V = 0 case, the preceding analysis reveals an issue at its core: the tension between the weak-noise and long-time limits, which do not commute [23]. Indeed, the Laplace parameter $s = H/4D$, and care is required with the standard Laplace transform notion of 'small $s\leftrightarrow$ large t' when D is small. A path involving multiple circuits of the interval will have a very large W, and one might be tempted to neglect it. However, it is required for the correct evaluation of the long-time limit—no matter how improbable a path may be made by the weakness of the noise, at sufficiently large times it will contribute. A partial explanation lies in the time-domain action, equation (15): the term −Ht acts to compensate a large-HW value (which $\sim 2\sqrt H$ at large H), increasingly so when t becomes large as well. This suggests that, in the H-domain, the important criteria for inclusion is that the path must satisfy the Euler–Lagrange equations for the stochastic action, rather than be a small-action path per se.

One way to approach simple diffusion equations with finite boundary conditions is to use the method of images, and there is a one-to-one correspondence between the bouncing paths considered here and image terms. The paths which bounce off x = a and b can be identified as images located just beyond the boundaries, and the series of higher-order paths are the images of the images. The path integral approach offers considerable scope for generalization to nonzero potentials, where the method of images does not obviously apply [29]. The path integral approach also offers a weak-noise approximation to the Green's function for the Smoluchowski equation, which could be applied numerically to more complicated geometries.

In this section, we apply the results above to the harmonic oscillator, with potential $V = \frac1{2}\alpha x^2$ and infinite boundary conditions $\rho(x,t)\to 0$ as $x\to\pm\infty$. We take the starting point x₀ = 0 for simplicity. The time-domain solution is the well-known Ornstein-Uhlenbeck (OU) density,

$$\begin{align} \rho_{\mathrm{OU}}\left(x,t\right) = \sqrt{\frac{\alpha}{2\pi D\left(1-{\mathrm{e}}^{-2\alpha t}\right)}}\exp\left(\frac{-\alpha x^2}{2D\left(1-{\mathrm{e}}^{-2\alpha t}\right)}\right). \end{align} \tag{ 28 }$$

Applying equation (17), the integrals can be performed explicitly, giving

$$\begin{align} W\left(x;H\right) = \alpha x^2+\frac1{\alpha}\left(\alpha x\,\sqrt{H+\alpha^2 x^2} - H\log\left( \frac{-\alpha x + \sqrt{H+\alpha^2 x^2}}{\sqrt H}\right)\right), \end{align} \tag{ 29 }$$

and hence

$$\begin{eqnarray} \bar\rho\left(x;H\right) & = & {\mathcal J}\exp\left(-\frac{W\left(x,H\right)}{4D}\right)\left(H\left(H+\alpha^2 x^2\right)\right)^{-1/4}\nonumber\\ \rho\left(x,t\right) & = & \frac1{2\pi i}\int_{\mathrm{Br}}{\mathcal J}\exp\left(+\frac{Ht}{4D}\right)\exp\left(-\frac{W\left(x,H\right)}{4D}\right)\left(H\left(H+\alpha^2 x^2\right)\right)^{-1/4} \frac{{\mathrm{d}} H}{4D} \end{eqnarray} \tag{ 30 }$$

for $\bar\rho$ and ρ, where the Laplace inversion is performed along the Bromwich contour. The exponent is not quadratic in H, so the method of steepest descents will not give the integral exactly. Moreover, the expression $\bar\rho$ does not solve the transformed problem exactly either (appendix B), so even if we could perform the inversion integral exactly, we would not expect to recover the OU solution. Proceeding optimistically with steepest descents, the stationary point of the exponent is at $H_*$, defined by equation (10)

$$\begin{align} t = \frac{\partial W}{\partial H}\bigg|_{H = H_*} = \int_0^x\frac{{\mathrm{d}} y}{\sqrt{H_* + \alpha^2y^2}} = \frac1{2\alpha}\log\left[ \frac{\alpha x+\sqrt{H_* +\alpha^2 x^2}}{-\alpha x+\sqrt{H_* +\alpha^2 x^2}}\right], \end{align} \tag{ 31 }$$

which can be inverted to give

$$\begin{align} H_* & = \alpha^2x^2\,{\mathrm{cosech}}^2\alpha t\in \mathbb{R}^+. \end{align} \tag{ 32 }$$

$$\begin{align} W\left(H_*\right) & = \alpha x^2\left(1+{\mathrm{coth}}\,\alpha t\right) + \alpha^2 x^2\,t\,{\mathrm{cosech}}^2\alpha t = \alpha x^2\left(1+{\mathrm{coth}}\,\alpha t\right) +H_*t, \end{align} \tag{ 33 }$$

so

$$\begin{align} S\left(x,t\right) = \alpha x^2\left(1+{\mathrm{coth}}\,\alpha t\right) = 2\alpha x^2/\left(1-{\mathrm{e}}^{-2\alpha t}\right). \end{align} \tag{ 34 }$$

The contour is then the line $H = H_*+iy$, $y\in \mathbb{R}$, and steepest descents gives

$$\begin{align} \rho\left(x,t\right) = \frac{{\mathcal J}\exp\left(-S/4D \right)}{\left(H_*\left(H_*+\alpha^2 x^2\right)\right)^{1/4}}\sqrt{\frac1{8\pi D|W^{^{\prime\prime}}\left(H_*\right)|}}. \end{align} \tag{ 35 }$$

The integral for W'' (see equation (14)) can also be done explicitly (see appendix A), and after some simple manipulations and the use of ${\mathcal J} = \exp(\alpha t/2)$ ¹ , this reduces exactly to the OU form equation (28). So, an approximate inversion of an approximate solution in the Laplace domain results in the exact time-domain solution. Of course, it would be too much to hope that the exactness would extend to more complicated potentials, but this does illustrate the surprising efficacy of the semiclassical approximation. One possibility is that this surprise stems from the underlying simplicity of the harmonic oscillator decoupling the fluctuations from the mechanical force, as observed in section 3.2 of [17].

In this paper we have revisited the path integral formulation of stochastic processes. A closer look at the correspondence between the dominant stochastic paths and conservative Hamiltonian mechanics in an effective potential allows the conserved 'energy' H (which has dimensions of power) to be identified as the Laplace parameter of the diffusive dynamics. $H\to 0,\,t\to\infty$ is interpreted as the equilibrium limit. H = 0 paths satisfy $\dot x = \mp V^{^{\prime}}$, and correspond to either noiseless relaxation to a local potential minimum, or 'optimal' (probability-maximsing) hill-climbing respectively. As might be expected from comparing the Schrödinger equation with the Smoluchowski equation, the Laplace transform of a diffusive process plays an analogous role to the Fourier transform in quantum mechanics. The latter is governed by oscillatory wave-like dynamics, whereas the stochastic case has decaying relaxational modes. Working in the Laplace domain facilitates summing over trajectories involving multiple circuits of some interval, because the series of contributing terms to the transition probability density is geometric. Indeed, including an infinite series of such turning paths is essential for correctly evaluating the long-time, equilibrium limit. Finally, we discussed a curious error cancellation in the harmonic oscillator, where the approximate Laplace inversion of an approximate solution in the Laplace domain led to the exact time domain solution.

S P F acknowledges many useful discussions with Profs A Archer and A McKane. This work was supported by the UK EPSRC, Grant Number EP/R005974/1.

No new data were created or analysed in this study.

In this appendix, we derive the explicit form of det$\,\hat L$, and discuss its long-time limit in the case of the harmonic oscillator. A convenient form for the determinant is given by (see e.g. [12, 13])

$$\begin{align} {\mathrm{det}}\,\hat L = \left|\text{det}\,\frac{\partial x_{*i}\left(t\right)}{\partial\dot x_{*j}\left(0\right)}\right| = \left|\frac{\partial x_*\left(t\right)}{\partial\dot x_*\left(0\right)}\right|\;\;\text{in 1D}. \end{align} \tag{ A.1 }$$

It is the derivative of the final position with respect to the initial momentum, and as such can be thought of as a 'density of paths'. On an extremal path $x_*$, $H = \dot x_*^2 - V^{^{\prime}} {}^2$ is constant. Dropping the ${}_*$ subscript and labelling $x(0) = x_0,x(t) = x_1$, evaluating H at t = 0 gives

$$\begin{align} \frac{\partial}{\partial \dot x_0} = \frac{\partial H}{\partial \dot x_0} \frac{\partial}{\partial H} = 2\dot x_0 \frac{\partial}{\partial H} = 2\sqrt{H+V^{^{\prime}}\left(x_0\right)^2} \frac{\partial}{\partial H}; \end{align} \tag{ A.2 }$$

(initial momentum and position are independent). Varying $H\to H+\delta H$, while holding t constant $(\delta t = 0)$,

$$\begin{eqnarray} t+\delta t & = & \int_{x_0}^{x_1+\delta x_1}\frac{{\mathrm{d}} x}{\sqrt{H+\delta H + V^{^{\prime}} {}^2}}\nonumber\\ & = & \int_{x_0}^{x_1+\delta x_1}\frac{{\mathrm{d}} x}{\sqrt{H + V^{^{\prime}} {}^2}} \left(1 -\frac1{2}\frac{\delta H}{H+V^{^{\prime}} {}^2} + {\ldots}\right)\;\;\;\nonumber\\ & \approx & t + \delta x_1 \frac1{\sqrt{H + V^{^{\prime}}\left(x_1\right)^2}} -\frac{\delta H}{2}\int_{x_0}^{x_1}\frac{{\mathrm{d}} x}{\left(H+V^{^{\prime}} {}^2\right)^{3/2}} , \end{eqnarray} \tag{ A.3 }$$

so to keep t constant, we require

$$\begin{align} \frac{\delta x_1}{\delta H} = \sqrt{H + V^{^{\prime}}\left(x_1\right)^2} \cdot \frac1{2}\int_{x_0}^{x_1}\frac{{\mathrm{d}} x}{\left(H+V^{^{\prime}} {}^2\right)^{3/2}}, \end{align} \tag{ A.4 }$$

and hence

$$\begin{align} \frac{\partial x_1}{\partial \dot x_0} = \sqrt{H+V^{^{\prime}}\left(x_0\right)^2} \sqrt{H + V^{^{\prime}}\left(x_1\right)^2} \int_{x_0}^{x_1}\frac{{\mathrm{d}} x}{\left(H+V^{^{\prime}} {}^2\right)^{3/2}}. \end{align} \tag{ A.5 }$$

If V has a critical point on the path, then the integral in the above diverges as H → 0. This is a manifestation of the well-known zero mode problem encountered when taking the long-time limit of fluctuation determinants (see e.g [28]). This is cancelled by the Jacobian, at least the in case of the harmonic oscillator. Taking $V = \frac1{2}\alpha x^2$ and x₀ = 0 for simplicity,

$$\begin{align} \int_{0}^{x_1}\frac{{\mathrm{d}} x}{\left(H+V^{^{\prime}} {}^2\right)^{3/2}} = \frac{x_1}{H\sqrt{H+\alpha^2 x_1^2}}\to\infty\;\text{as}\; H\to 0, \end{align} \tag{ A.6 }$$

and equation (32) gives

$$\begin{align} H = \alpha^2x_1^2\,{\mathrm{cosech}}^2\alpha t;\;\;\sqrt{H+\alpha^2 x_1^2} = \alpha x_1\,{\mathrm{coth}}\,\alpha t. \end{align} \tag{ A.7 }$$

This results in

$$\begin{align} \frac{\partial x_1}{\partial \dot x_0} = \frac{\sinh\,\alpha t}{\alpha};\;\; \rho = \sqrt{\frac{\alpha}{\sinh\,\alpha t}}\,{\mathcal J}\exp{-\left(S/4D\right)}. \end{align} \tag{ A.8 }$$

Were it not for ${\mathcal J}$, this would → 0 as $t\to\infty$, but the initially somewhat worrying-looking ${\mathcal J} = \exp\left( +\frac1{2}\alpha t\right)$ means the correctly-normalised OU density equation (28) is recovered, with the long-time limit

$$\begin{align} \rho_{\mathrm{eq}}\left(x\right) = \sqrt{\frac{\alpha}{2\pi D}}\exp\left(-\frac{\alpha x^2}{2D}\right). \end{align} \tag{ A.9 }$$

The Laplace-transformed Smoluchowksi equation is

$$\begin{align} s\bar\rho\left(x;s\right) - \delta\left(x-x_0\right) = \frac{{\mathrm{d}}}{{\mathrm{d}} x}\left(\bar\rho V^{^{\prime}}+D\frac{{\mathrm{d}} \bar\rho}{{\mathrm{d}} x}\right). \end{align} \tag{ B.1 }$$

Inserting the WKB ansätz

$$\begin{align} \bar\rho = \exp\left( -\frac1{D}\sum_{n = 0}^{\infty} D^nf_n\left(x\right)\right) = A\left(x\right)\exp\left(-\frac{W\left(x\right)}{4D}\right) + \text{h.o.t.} , \end{align} \tag{ B.2 }$$

and keeping the first two terms only, leads to

$$\begin{align} As \approx V^{^{\prime}}A^{^{\prime}}-\frac{AV^{^{\prime}}W^{^{\prime}}}{4D}+AV^{^{\prime\prime}}+DA^{^{\prime\prime}}-\frac{AW^{^{\prime\prime}}}{4} -\frac{A^{^{\prime}}W^{^{\prime}}}{2}+\frac{AW^{^{\prime}} {}^2}{16D} . \end{align} \tag{ B.3 }$$

Identifying $s = H/4D$, and thus including it in the dominant balance, immediately gives

$$\begin{align} W^{^{\prime}} = 2V^{^{\prime}}\pm 2\sqrt{H+V^{^{\prime}} {}^2};\; W = 2V\pm2\int\sqrt{H+V^{^{\prime}} {}^2}\,{\mathrm{d}} x \end{align} \tag{ B.4 }$$

exactly as with the path integral, with ± corresponding to right and left-moving paths respectively. If we had excluded the Laplace parameter from the dominant balance, we would have found $W = 4V$, which is only correct in the long-time, time-independent equilibrium limit. The ${\mathcal O} (1)$ terms give

$$\begin{eqnarray} 0 & = & A^{^{\prime}}\left(V^{^{\prime}}-\frac{W^{^{\prime}}}{2}\right) + A\left(V^{^{\prime\prime}}-\frac{W^{^{\prime\prime}}}{4}\right); \nonumber\\ \log A & = & -\frac1{4}\log\left(H+V^{^{\prime}} {}^2\right) \pm\frac1{2}\int^x V^{^{\prime\prime}}\frac{{\mathrm{d}} x}{\sqrt{H+V^{^{\prime}} {}^2}}\nonumber\\ A & = & \left(H+V^{^{\prime}} {}^2\right)^{-1/4} \exp\left(+\frac1{2}\int_0^tV^{^{\prime\prime}}\left(x\left(t\right)\right){\mathrm{d}} t\right), \end{eqnarray} \tag{ B.5 }$$

where the last line follows using $\dot x = \pm\sqrt{H+V^{^{\prime}} {}^2}$. (Within the steepest descents approximation, since ${\mathcal J}$ is independent of the small noise D, it is not involved in the integral, and can be thought of simply as ${\mathcal J}(t(H_*))$ in the Laplace domain). So, the 1/4 power prefactor term and the Jacobian emerge from a truncated WKB expansion, as expected by analogy with the quantum-mechanical case. The path integral formula gives the correct normalizing factor, and prescription for left and right-moving paths, whereas the full WKB solution would have to be constructed by explicitly matching solutions,and imposing the appropriate jump condition on $\bar\rho^{^{\prime}}$ at $x = x_0$. However, even for the harmonic potential $V = \alpha x^2/2$, this is not exact at this order: the term DA'' is neglected. By approximately inverting the transform using steepest descents, the exact time domain (OU) solution is recovered.