{'fraction_non_alphanumeric': 0.042746909459856354, 'fraction_numerical': 0.014818572611005793, 'mean_word_length': 4.260393258426967, 'pattern_counts': {'":': 0, '<': 14, '<?xml version=': 0, '>': 12, 'https://': 1, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 9, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
The computational cost of blocking for sampling discretely observed diffusions

Many approaches for conducting Bayesian inference on discretely observed diffusions involve imputing diffusion bridges between observations. This can be computationally challenging in settings in which the temporal horizon between subsequent observations is large, due to the poor scaling of algorithms for simulating bridges as observation distance increases. It is common in practical settings to use a blocking scheme, in which the path is split into a (user-specified) number of overlapping segments and a Gibbs sampler is employed to update segments in turn. Substituting the independent simulation of diffusion bridges for one obtained using blocking introduces an inherent trade-off: we are now imputing shorter bridges at the cost of introducing a dependency between subsequent iterations of the bridge sampler. This is further complicated by the fact that there are a number of possible ways to implement the blocking scheme, each of which introduces a different dependency structure between iterations. Although blocking schemes have had considerable empirical success in practice, there has been no analysis of this trade-off nor guidance to practitioners on the particular specifications that should be used to obtain a computationally efficient implementation. In this article we conduct this analysis and demonstrate that the expected computational cost of a blocked path-space rejection sampler applied to Brownian bridges scales asymptotically at a cubic rate with respect to the observation distance and that this rate is linear in the case of the Ornstein-Uhlenbeck process. Numerical experiments suggest applicability both of the results of our paper and of the guidance we provide beyond the class of linear diffusions considered.


Introduction
Diffusions have been widely applied to model continuous-time phenomena of interest, including molecular dynamics (Boys et al. 2008), neuroscience (Lansky & Ditlevsen 2008), and finance (Karatzas & Shreve 1998). In general, a diffusion on R d is a Markov process X defined to be the solution, with law we will denote by P, to a stochastic differential equation of the following form: where b ∶ R d → R d and σ ∶ R d → R d×d ′ denote the drift and volatility coefficient respectively, and W is a standard d ′ -dimensional Brownian motion. Throughout we assume standard regularity conditions hold which ensure the existence of a unique, global, weak solution to (1) (see for instance Øksendal 2003).
In practice we will typically only have access to discrete observations of (1), and so for practitioners the statistical problem of interest is to use these observations to draw inference on the parameters of b and σ of (1). A common Bayesian strategy is to augment the parameter space with the space describing the complete underlying diffusion trajectory. A Markov Chain Monte Carlo algorithm can then explore this augmented space by alternating between updates of the parameters and updates of the unobserved sample path connecting observations (sampling of diffusion bridges) (Roberts & Stramer 2001). As a consequence a considerable and methodologically diverse literature has been developed concerned with simulating diffusion bridges (the law of (1) conditioned to terminate at the subsequent observation (for instance, X T = x T ), which we denote by P (T,x0,x T ) or generically by P ⋆ ), including Beskos & Roberts (2005), Bladt et al. (2014), Delyon & Hu (2006), Durham & Gallant (2002), Golightly & Wilkinson (2008), Hairer et al. (2011), Roberts & Stramer (2001), Schauer et al. (2017).
One of the common difficulties with Markov Chain Monte Carlo strategies is sampling diffusion bridges between distant observations; the duration of the bridge, which we denote by T , is large. This setting naturally arises when the underlying diffusion (1) is sparsely observed (or high-dimensional), for instance in shape analysis applications (Arnaudon et al. 2020), or in the case of diffusions on graphs (Freidlin & Wentzell 1993). The problem here is that methodologies for sampling diffusion bridges scale poorly with T , and many of the most widely used approaches have exponential computational cost in T . Consequently, addressing the poor scaling in T has drawn considerable interest. One popular approach is the blocking scheme introduced by Shephard & Pitt (1997), which has been employed in a number of practical problems with strong empirical evidence of its efficacy (Chib et al. 2004, Golightly & Wilkinson 2008, Kalogeropoulos 2007, Kalogeropoulos et al. 2010, van der Meulen & Schauer 2018, Stramer & Roberts 2007. Blocking is a conceptually simple idea in which the time domain of the diffusion bridge is overlaid with a set of temporal anchors (0 =∶ k 0 < k 1 < ⋅ ⋅ ⋅ < k m < k m+1 ∶= T ), and the values of the bridge are taken for some initialisation trajectory at those points (which are known as knots, and for which we will denote X i ∶= X ki to simplify notation).
Simulation from P (T,x0,x T ) is then achieved by constructing a Gibbs sampler which alternates between updating knots and updating the segments of the trajectory conditional on the knots, a number of times. For instance, we could begin by simulating from the conditional law P (k2−k0,X0,X2) (updating the trajectory between [t 0 , t 2 ] which includes the knot at X 1 , conditional on the knots at X 0 and X 2 ), and then P (k3−k1,X1,X3) (updating the trajectory between [t 1 , t 3 ] and containing the knot X 2 , conditional on the knots at X 1 and X 3 ), and so on sweeping across all anchor points. This sweep would then be iterated a number of times to reduce the dependency between the resulting bridge and the initial (or previous) trajectory. In this article we consider the three canonical blocking schemes of Roberts & Sahu (1997) with equidistant anchors: the checkerboard scheme, in which the odd and even indexed knots are alternatively updated; the lexicographic scheme, in which the knots are updated in temporal order; and the random scheme, in which at each step a random knot is updated. We more formally introduce blocking and define these schemes in Section 2.
From a computational perspective, blocking substitutes the expensive simulation of a (single independent) draw from P (T,x0,x T ) , with the cost of simulating repeated sweeps of the m + 1 shorter (and computationally more efficient) bridges for each segment given by the temporal anchors. Any analysis of this trade-off needs to take into account the serial correlation induced by the blocking strategy.
Despite widespread adoption of blocking in practice to mitigate the computational cost of simulating diffusion bridges (as indicated above), there is little theoretical support for its efficacy. Furthermore, there is little concrete guidance on how to implement, and then appropriately tune (selecting for instance the number and locations of the anchor points), a blocking scheme.
In this article we provide general guidance for implementing blocking schemes by addressing these practical considerations for particular classes of diffusion process. We analyse the computational cost of several rejection sampling algorithms for bridges as a function of block size and bridge duration. In all cases we consider a fixed regular spacing of m anchor points as m, T → ∞, in contrast with the study of the 'in-fill' asymptotic of Roberts & Stramer (2001) in which T is fixed and m → ∞. We analyse the expected cost of a single iteration of various algorithms, and then to capture the trade-off described above we consider the cost of the algorithm which comprises both the cost of one iteration, and the total number of iterations required to obtain an 'independent' sample. We give a more formal description of what we mean by achieving independence below, in terms of the relaxation time of the underlying Markov chain.
In this article we work under the assumption that the underlying measure is a Gaussian diffusion (i.e. P is the law of a scaled Brownian motion or the law of the Ornstein-Uhlenbeck process). Under this simplification the Gibbs step for updating the bridge segments can be implemented without error, i.e. without discretising time, for example by means of a rejection sampler directly on the path-space of the diffusion (see Appendix A for full details). In this setting we prove that Theorem 1 below holds, as the culmination of the results in Section 3. We gather all proofs in the appendices. Theorem 1. Suppose P ⋆ is the conditional law of a Gaussian diffusion which is sampled by rejection on path-space and using a checkerboard or lexicographic or random blocking scheme. Suppose the m anchors are spaced equidistantly such that m = c 1 T (for some constant c 1 > 0). Then the expected computational cost of the blocked rejection sampler, C blocking (T ), satisfies: whenever P denotes the law of a scaled Brownian motion and whenever P denotes the law of the Ornstein-Uhlenbeck process.
Remark. Note that in the case of a Brownian bridge there is long range dependency in the path, in the sense that the correlation between X s and X t is non-negligible even for 0 ≪ s ≪ t ≪ T . On the other hand, for the Ornstein-Uhlenbeck process its ergodicity breaks this dependency. For an Ornstein-Uhlenbeck process whose drift is of the form b(X t ) = −θX t and for T ≫ 0, there is a phase transition in its behaviour as θ → 0 in that the computational cost of a blocked rejection sampler for Brownian motion is not recovered. Recall that in this paper we are working under the assumption that the underlying measure is a Gaussian diffusion, but in most practical settings the target law will be more complicated. In such settings it would be typical to use a Gaussian diffusion as a proposal law for the non-Gaussian target law. In principle Theorem 1 would suggest that an Ornstein-Uhlenbeck process proposal for a stationary target law would be advantageous over a Brownian bridge proposal, although in practice this predicted computational saving would depend on how closely the target process matched the invariant distribution of the Ornstein-Uhlenbeck process.
Theorem 1 contrasts sharply with the case without blocking. We show later in Proposition 1 that, for a d-dimensional Brownian bridge proposal in the absence of blocking, the cost is exponential in T . Although what we prove in Theorem 1 addresses a somewhat idealised setting, the requirement m = c 1 T acts as a concrete guide for choosing the number of blocks. Furthermore, our empirical results in Section 4 indicate that the guidance we establish can be more broadly useful beyond the class of linear diffusions. Thus we demonstrate that blocking can lead to significantly improved computational efficiency when conducting inference for discretely observed diffusions.

Blocking
In this section we provide a systematic definition of blocking for sampling a diffusion path. Define a set of anchors spread across the time domain: 0 < k 1 < ⋅ ⋅ ⋅ < k m < T and knots as the values of the path taken at the anchors: Each anchor is now assigned to one of k disjoint subsets A i , i = 1, . . . , k, each comprising m i anchors: In particular, A i = {r i1 , . . . , r imi } is the set of anchors associated (uniquely) to the ith block. This allows us to group the knots by associating them with the corresponding subsets of anchors: . . , k}.
For convenience of notation we let K −i (resp. A −i ) denote all knots (resp. anchors) that do not belong to K i (resp. A i ): and assign labels to an ordered collection of all anchors in A −i , plus the end-points: , to be only those intervals between the end-points or anchors in A −i , which contain at least one anchor belonging to A i .
The path segments X Bi , obtained through restricting X to B i , are termed blocks. Finally, in the case k = 2 we say that A 1 and A 2 are interlaced if whenever a, c ∈ A i , with a < c, then there exists b ∈ A (i mod 2)+1 s.t. a < b < c, i = 1, 2.
A sampler for a path equipped with a blocking technique is a Gibbs sampler that updates the full path only one block at a time by drawing from the conditional laws P ⋆ Bi (⋅ K −i )-i.e. the target laws restricted to blocks B i and conditioned on the knots in K −i . For simplicity we refer to this technique as a blocked sampler in the remainder of the text, and present general pseudo-code for it in Algorithm 1.
Algorithm 1: Blocked sampler on path-space There are a number of ways we can update the blocks, and in this article we consider the three canonical blocking schemes of Roberts & Sahu (1997). In particular, we refer to a single, full Gibbs sweep of Algorithm 1 (the inner for-loop) as a: Definition 1. Checkerboard blocking update scheme if k = 2, A 1 and A 2 are interlaced, and q(i, j) ∶= 1 {i} (j).
The above are not exhaustive, but characterise the most widely used, and are tractable enough for analysis. We further simplify various computations by assuming the anchors are equidistant, and defer discussion of this assumption and its relaxation to Section 5. Assumption 1. The anchors are placed on an equidistant grid: As mentioned in the Introduction, we will study the asymptotic regime in which δ m,T is fixed as m, T → ∞.
3 Computational analysis

Cost of a single sweep
We begin by quantifying the computational cost of a rejection sampling algorithm for diffusion bridges in the absence of blocking. The setting here, as given in further detail in Appendix A, is one in which the target law is absolutely continuous with respect to a d-dimensional Brownian bridge proposal path. Proposition 1. Under Assumptions 3-9 enumerated in Appendix A, the expected computational cost as a function of T of obtaining a single draw with a path-space rejection sampling algorithm, denoted by C rej (T ), is given by where c 2 > 0 is some constant independent of T , and the function f ∶ R + → R is continuous and such that f (T ) ∼ T −d 2 as T → ∞. In particular, for large enough T there is a constant c 3 > 0 such that: Remark. Note that Proposition 1 does not stipulate in what way the constant c 3 might vary with dimension. Without further structure it is impossible to characterise this behaviour. However it is highly likely that c 3 will increase at least linearly with dimension (as for instance would be the case for diffusions consisting of d independent components).
Remark. If P is the law of a drifted Brownian motion, then Proposition 1 cannot be applied directly, because Assumption 9 does not hold. However, for this case an easy calculation shows that the acceptance probability of a rejection sampler with Brownian bridge proposals is equal to 1, implying (under Assumption 8) that C rej (T ) is proportional to T . Now considering a single sweep of the blocking schemes introduced in Section 2, note that we have substituted sampling a single diffusion bridge (of length T ) with sampling a number of diffusion bridges of shorter time horizon, 2δ m,T (for example, to ensure the point X 2δ m,T is updated one could sample a new bridge of length 2δ m,T connecting X δ m,T with X 3δ m,T ). By application of Proposition 1, the expected computational cost of simulating each of these shorter bridges is therefore C rej (2δ m,T ), and hence the expected cost of a single Gibbs sweep is: Equation (5) holds for all m and T and follows from (4); however, as the behaviour of f (t) for small t is not immediately transparent, to learn something about C sweep (T, m) when δ m,T is small, we may use the fact that the acceptance probability of the rejection sampler approaches 1 as the bridge duration decreases to 0. This fact implies that for small enough t, C rej (t) ∼ c 5 t and thus if δ m,T < c 4 as T → ∞ for some constant c 4 then for some c 5 > 0. For instance, upon setting m = ⌊T ⌋, the cost in (5) becomes O(T ) as T → ∞. Contrast this with (4) to see that the relative gain in efficiency, C rej (T ) C sweep (T, m) grows exponentially in T and suggests that blocking is to be preferred for large enough T . However, this ignores the costs associated with mixing; we address this in the next subsection.

Cost of multiple sweeps
Direct comparison of the exponential cost C rej (T ) of direct rejection sampling (as given by Proposition 1), with the linear cost C sweep (T, m) of a single sweep of a blocking scheme (as given by (6)), does not capture the remnant dependency structure introduced by the blocking scheme. In addition we need to consider the number of sweeps required to render this dependency negligible. In order to do that we first introduce the following notion Definition 4. (Roberts & Sahu 1997) The [L 2 -]convergence rate ρ of a Markov chain {X (n) ; n = 1, . . . , N } with the transition kernel P and an invariant density π is defined as the minimum number for which for all square π-integrable functions f , and for all r > ρ We can now capture the cost of reducing the dependency on the past by considering the relaxation time, denoted T = T (T, m), and defined as: It represents the time required by the underlying Markov chain to output a draw from its stationary distribution (Levin & Peres 2017). This makes it possible to compare C rej (T ) with the expected computational cost of the blocked rejection sampler as follows: We will later consider the most appropriate choice of blocking scheme, and how to optimise m.
Instead of analysing the chain targeting the law P ⋆ it is sufficient to consider a related chain that targets the marginal law of the vector G ∶= (X k1 , . . . , X km ) (X 0 , X T ) = K (X 0 , X T ), which we denote by G. To see this, notice that conditionally on the knots K being distributed according to G, a path X returned after a single Gibbs sweep of a blocking scheme is distributed exactly according to P ⋆ . The object of interest becomes a Markov chain with a transition kernel P denoting a single Gibbs sweep, and with stationary distribution Throughout, we additionally assume that the following condition holds, which makes the subsequent required calculations tractable.
Assumption 2. The target law P ⋆ is such that G is a Gaussian process.
We discuss this key technical assumption in Section 4, where we note that the established results seem to hold empirically more broadly.
Under Assumption 2 and using either the lexicographic or checkerboard updating scheme, a single Gibbs step (i.e. an update G B I ∼ P ⋆ B I ∩{X k 1 ,...,X km } (⋅ K −I )) has a tractable, Gaussian transition density, and thus so does the entire Gibbs sweep G (n) ↦ G (n+1) with mean and covariance As a consequence it is possible to explicitly characterise the transition kernel P , as follows.
Lemma 1. Under the lexicographic and checkerboard updating schemes, the n-step transition kernel P n of the Markov chain {G (l) ; l = 0, . . . } is Gaussian, with mean and covariance matrix given respectively by: with B ∈ R m×m and b ∈ R m .
Under the lexicographic or the checkerboard updating schemes {G (l) ; l = 0, . . . } is an AR(1) process, and so the spectral radius ρ spec (B) of the matrix B must satisfy ρ spec (B) < 1 for the process to converge, and equals the L 2 -convergence rate (Amit 1991). This connection extends to the random updating scheme. In the following lemma we derive the spectral radius of each blocking scheme as a function of m and T , which aids in optimising their parameterisation and analysing their scaling. We denote by Λ ∶= Σ −1 the precision matrix of G and define A ∶= I − diag{Λ −1 11 , . . . , Λ −1 mm }Λ. Lemma 2. (Roberts & Sahu 1997) Under the checkerboard and lexicographic updating schemes, the spectral radius of the matrix B and the L 2 -convergence rate of a blocked rejection sampler coincide. More explicitly, under the checkerboard, lexicographic, and random updating schemes respectively the L 2 -convergence rates (ρ check , ρ lex , and ρ rand resp.) are equal to: where λ max (A) denotes the maximum eigenvalue of the matrix A and where we write B check (resp. B lex ) to denote a matrix B corresponding to the checkerboard (resp. lexicographic) updating scheme.
λ max (A) can be found more explicitly by exploiting the close connection between the precision matrix Λ and the matrix of partial correlations (given precisely in (B.2), in Appendix B).
with c(δ m,T ) ∶= Corr(X δ , X 2δ X 0 , X 3δ ). In particular: The form of c(δ m,T ) will, in general, depend on the type of a Gaussian process that is being considered. In the following corollaries we present more explicit versions of the statements from Theorem 2 for the two choices of P: scaled Brownian motion σW , with σ > 0; and, the Ornstein-Uhlenbeck process. Without loss of generality we centre the latter at 0: Corollary 1. If P is the law of a scaled Brownian motion σW , σ > 0, then: In particular, independently of T , as m → ∞ Corollary 2. If P is the law of the Ornstein-Uhlenbeck process (12), then: In particular, when δ m,T = δ is set to a constant, as m, T → ∞: Remark. Results obtained by Pitt & Shephard (1999), who studied the discrete-time first-order autoregressive process α t = φα t−1 + η t , η t ∼ N (0, σ 2 ), observed with Gaussian noise, are closely related to Corollaries 1 and 2. In the context of their model, where c(δ m,T ) = φ (1 + φ 2 ), they derive the expression (11) for ρ m,T as well as bounds on ρ m,T which exhibit the same asymptotic behaviour as in Corollary 2.
We can now combine the above results with (7) to find the relaxation time: Theorem 3. Suppose we use one of the checkerboard, lexicographic, and random updating schemes. If P is the law of a scaled Brownian motion σW , then we have: If P is the law of the Ornstein-Uhlenbeck process in (12), and additionally the sequence T (m) is chosen so that m = c 1 T for some constant c 1 > 0, then we have: Remark. Note that if P is the law of the Ornstein-Uhlenbeck process in (12) then Theorem 3 holds for only the sequence T (m) where m = c 1 T , but in the case of scaled Brownian motion there is no such constraint; see Remark 1.
Remark. From the proof of Theorem 3, one can show that for the Ornstein-Uhlenbeck process (12): This provides insight into the influence of θ and δ on mixing.
We can minimize the cost of blocking C blocking (T, m) over the remaining parameter, m, using Theorem 3 and (8). This leads to Theorem 1, which is the main result of this paper (as presented in Section 1, with accompanying proof in Appendix B).

Numerical experiments
Consider a target process defined to be the solution of the following stochastic differential equation (with law P): This diffusion exhibits highly multimodal behaviour, and so in practice it is challenging to simulate trajectories of P (and in particular the conditioned bridge law P ⋆ over large time horizons). It is possible to simulate trajectories exactly by means of path-space rejection sampling (as detailed in Appendix A). However, X is not a Gaussian process (it violates Assumption 2), and so Theorem 1 does not hold in a rigorous sense. As such (13) makes an interesting case to investigate the practical limitations of Theorem 1. Because it is not an ergodic diffusion, out of the two theoretical results from Theorem 1 the ones for the Brownian motion are expected to be more relevant. As we show below, the empirical results would suggest the theory holds more broadly.
We consider six problems (increasing in difficulty) of simulating paths according to the laws P (T,x0,x T ) , with parameters P (0.2,0,0.1) , P (0.4,0,0.85) , P (0.5,0,0.85) , P (1,0,0.95) , P (2,0,2.5) , P (4,0,4.85) . The values of the end-points were chosen by fixing T , simulating multiple paths according to (13) and picking x T to be some point in the vicinity of the (largest) mode as these are the bridges we will most commonly be interested in. For T = 0.2, the plotted paths resemble Brownian bridges, but as T increases the non-linear dynamics become pronounced: the diffusion is effectively attracted to a ladder of values and it is repelled at the intermediate points, leading to multimodal behaviour of the trajectories. Drawing paths from the last three laws using path-space rejection sampling but without blocking (an unmodified rejection sampler) is computationally infeasible.
For each of the six examples we ran a blocked rejection sampler with checkerboard updating scheme for 10 5 iterations and with various numbers of knots. For the first three problems we also employed an unmodified rejection sampler. We recorded the time required to sample a single path (which for a blocked rejection sampler is counted as one execution of the inner for-loop of Algorithm 1) and plotted it in Figure 1 against the number of used knots. Code sufficient for reproducing these results can be found at https://github.com/mmider/blocking.
For T = 0.2 the unmodified rejection sampler clearly outperforms any blocking scheme. This is unsurprising as paths under P (0.2,0,0.1) closely resemble Brownian bridges (and indeed every diffusion behaves as a drifted Brownian motion on a small-enough time-scale). However, as T increases, this pattern changes and blocking reduces the cost of obtaining any single sample path. In particular, notice a steep, exponential reduction in cost that is especially pronounced for  (and ESS was computed according to Gelman et al. (2013, Section 11.5)) against the (half-) length of blocks (i.e. δ m,T ). As defined in Figure 2, taESS is approximately equal to a number of independent samples that can be drawn in one second. Clearly, the larger taESS is the more efficient the algorithm is. First, for any experiment we expect there to be a point for which increasing the number of knots any further will only lead to a decrease in taESS-this corresponds to all costs being dominated by the cost due to a slowdown in mixing and it is clearly illustrated by sharp dips of curves on the left side of Figure 2. Second, for examples for which the target law is sufficiently different from the law of Brownian bridges we expect that some level of blocking will improve the overall computational cost. This is also confirmed by the declines of taESS curves toward the right side of Figure 2. We note that under the most difficult sampling regimes it was impractical to run the algorithm with even fewer blocks due to excessive execution times-had the examples been run and the curves continued, the decline in performance would have been even starker. Additionally, Figure 2 is suggestive of there being an optimal value of δ m,T (somewhere around δ m,T ≈ 0.1), that is almost independent of T and m and that yields the highest taESS in each experiment. This is consistent with the results of Section 3, where an optimal number of knots was found to be m = c 1 T for some c 1 > 0, which implies the claim about the dependence of the optimal δ m,T on T and m.
Finally, we verify the bound from (2) empirically. To this end, notice that taESS −1 is approximately equal to the amount of time needed to obtain a single independent sample. This is consistent with the characterisation of the computational cost of a blocked rejection sampler as given in (8). Theorem 1 asserts that this cost scales at a cubic rate in the duration of the bridge, so long as δ m,T is set to a constant when T → ∞. Consequently, taESS(T ) should be at most a cubic function of T and if plotted on a log-log scale, this would be equivalent to taESS(T ) tracking some line with slope 3. Figure 3 gives this precise plot, showing that the prediction (2) is indeed satisfied.

Discussion
In this article we have analysed and provided practical guidance for using blocking schemes when conducting Bayesian inference for discretely observed diffusions. We achieved this by studying the computational cost of diffusion bridge sampling algorithms. We have shown rigorously that the computational cost of rejection sampling on path-space (modified with blocking) targeting the law of scaled Brownian motion scales as O(T 3 ) as T → ∞, and as O(T ) in the case of the Ornstein-Uhlenbeck process, so long as the number of equidistant anchors is m = c 1 T (for some c 1 > 0). In Remark 1 we discussed the practicality of exploiting the computational saving achievable in the case of the Ornstein-Uhlenbeck process. Furthermore, using the example of a non-linear sine diffusion we provide empirical evidence which would suggest that the conclusions about Brownian motion hold also for non-ergodic diffusions outside of a restrictive class of Gaussian processes.
Our theory indicates that choosing too few knots results in the computational cost being dominated by the exponential cost for imputing diffusion bridges between successive knots (see Proposition 1). As such our guideline of choosing m = c 1 T , (for some c 1 > 0) is useful for ensuring the robustness of blocking schemes and a reasonable heuristic for practitioners. Note that although choosing too many knots is likely to be penalized less than choosing too few, choosing an excessive number of knots can negatively impact the mixing of the underlying chain.
Naturally, for more general target laws P it might be useful to consider using irregularly spaced anchors (and so relaxing Assumption 1). Heuristically, we may wish to place more knots in areas in which the proposal law does not approximate the target law well. Developing more general theory to support the use of an irregular spacing of anchors is likely to require more knowledge of the specific diffusion under study. Of course, from a methodological perspective this motivates future research looking at how to place knots by assessing proposal-target discrepancy, or developing adaptive schemes.
Finally, it is worth recalling that within the context of Bayesian inference for discretely observed diffusion processes, the full chain in this setting is a Gibbs sampler that alternates between updating the unknown parameters and imputing the unobserved path. Since the mixing time of the unobserved path influences the mixing time of the parameter chain, then in light of the work in this paper it may as a future extension also be possible to study the mixing behaviour of the parameter chain.

A Rejection sampling on path-space
In this article we have restricted our attention to the class of diffusion bridges which can be sampled by means of path-space rejection sampling. In particular, to sample from P ⋆ we sample trajectories from an accessible and absolutely continuous proposal law (denoted Q * ), and accept with probability proportional to the Radon-Nikodým derivative of P ⋆ to Q * (Beskos & Roberts 2005, Beskos et al. 2006, Pollock et al. 2016. To find an appropriate Q * we impose the following common assumption (Kloeden & Platen 2013, Section 4.4) Assumption 3. There exists η ∶ R d → R d such that ∇η = σ −1 .
Under Assumption 3 the process Y ∶= {η(X t ), t ∈ [0, T ]} satisfies the following stochastic differential equation, for a known, closed-form drift α. With unit volatility, the law of Y is now absolutely continuous with respect to Brownian motion, and so Brownian motion is a viable proposal law for sampling from P (and the law induced by the Brownian bridge is a viable proposal for P ⋆ ).
In order to avoid unnecessary inflation of notation we assume throughout the article that σ ≡ 1 in (1), so that α ≡ b, η becomes an identity map, and X ≡ Y . The general case of σ (that satisfies Assumption 3) follows without additional effort. Assumption 4. α is at least once continuously differentiable.
Assumption 5. There exists a potential function A ∶ R d → R such that ∇A = α.
Under assumptions 4-6 we have (Beskos & Roberts 2005, Section 3): It follows that sampling from P ⋆ can be accomplished using Algorithm 2. Note that computing the integral in (A.1) required for Algorithm 2 can be achieved either (i) approximately, by simulating a candidate Y ○ over a fine mesh and computing the integral in (A.1) numerically; or (ii) exactly, via an additional randomisation step that utilises a Poisson point process (Beskos & Roberts 2005, Beskos et al. 2006. We refer to these two methods as, respectively, approximate and exact path-space rejection samplers. If we contrast U ∼ Bern(p(Y )) with that ofŨ ∼ Bern(p(Y )) wherep is the (unbiased) estimator resulting from the additional randomisation step, then since E[p] = p we have U d To prove Proposition 1, which quantifies the computational cost of Algorithm 2, we impose the following natural assumptions on the cost of simulating a proposal trajectory. Assumption 7. The cost of generating any proposal sample X is independent of the value of p(X) as defined in (A.1).
Assumption 8. The cost c(X) of simulating X has expectation growing linearly in T: E [c(X)] = c 5 T , c 5 ∈ R + .
Assumptions 7 and 8 are always satisfied if rejection sampling is performed with the approximate method described above, so long as the mesh width is kept constant as T → ∞. For the exact method, Assumption 7 will in general be violated (for instance, if the number of simulated Poisson points is 0, then conditional on this information p(X) = 1 a.s.) but for T → ∞ it is a reasonable approximation. Assumption 8 is satisfied if φ is bounded.
We can now derive the cost of a single draw using a path-space rejection sampler as follows: Lemma 3. Under Assumptions 3-8: where c 6 > 0 is some constant independent of T , p T (x 0 , x T ) is the transition density under P for going from x 0 to x T over the interval [0, T ] and q T (x 0 , x T ) is the same transition density, but under the proposal law Q instead.
Proof. Denote by X (i) , i ∈ {1, 2, . . . }, independent samples from Q * and by c(X (i) ) the cost of sampling path Rejection sampling requires a geometrically distributed number of simulations (with a randomly distributed parameter at each trial), so its expected cost is where the measures with respect to which the expectations above are taken should be clear from the context (and include Brownian bridge measures, products of Brownian bridge measures, and any additional randomness needed to simulate events of probability p(X (i) )). We now have: To better understand the scaling with T of the ratio of transition densities under the laws P and Q in Lemma 3 we impose the following final assumption, which allows us to establish Lemmata 4 and 5 required for proving Proposition 1.
Assumption 9. The target diffusion is ergodic and defined on R d .
Lemma 4. Under Assumption 9, for we have that f (T ) ∼ T −d 2 as T → ∞ and d denotes the dimension of the process.
Proof. q and p are well-behaved densities which are bounded and bounded away from zero (see for example Rogers (1985)). This implies that f is continuous. As the target diffusion is ergodic, p T (x 0 , v) →p(v) as T → ∞, wherep is the stationary density of the diffusion law. On the other hand q T (x 0 , v) is just a Gaussian density with variance T d I, which for T → ∞ behaves as ∼ T −d 2 .
Proof. From Lemma 4 we have that the RHS of (A.4) is ∼ T d 2 e ΦT for T → ∞. As the LHS of (A.4) represents an expected probability, we must have Φ < 0 for this expression to take values in [0, 1].
We are now in a position to prove Proposition 1.
Proof of Proposition 1. This follows directly by combining Lemmata 3, 4 and 5.

B Proofs
Proof of Lemma 1. The chain {G (l) ; l = 0, . . . } coincides with the chains considered in Roberts & Sahu (1997). In particular, the 1-step transition kernel under lexicographic and checkerboard updating schemes is stated explicitly as Roberts & Sahu (1997, Lemma 1). We provide a proof for completeness.
Proof of Lemma 2. Since the chain {G (l) ; l = 0, . . . } coincides with the chains considered in Roberts & Sahu (1997), the statement of Roberts & Sahu (1997, Theorem 1) applies under checkerboard and lexicographic updating schemes: i.e. the L 2 convergence rates under the two regimes are given by ρ check = ρ spec (B lex ) and ρ lex = ρ spec (B lex ) respectively. Due to tridiagonal structure of the precision matrix Λ (which follows from the Markov property of the process G; see also a short explanation in the proof of Theorem 2 that leads up to (B.3)), Roberts & Sahu (1997, Corollary 3) implies that the two spectral radii coincide, i.e. ρ spec (B lex ) = ρ spec (B check ). By the same token, Roberts & Sahu (1997, Theorem 5) applies as well, yielding ρ spec (B check ) = λ 2 max (A). Finally, the L 2 convergence rate of the random updating scheme follows from Roberts & Sahu (1997, Theorem 2).
Proof of Theorem 2. The precision matrix Λ of any random vector G with non-degenerate covariance matrix can be related to a matrix of partial correlations via (Lauritzen 1996, p. 130): j,j] .

(B.2)
By the definition of G in (9), it is easy to see that Corr(G The eigenvalues of Toeplitz matrices may be found in closed form (Smith 1985, Kulkarni et al. 1999) and in particular, those of matrix A are given by −2c(δ m,T ) cos πl m + 1 , l = 1, . . . , m.
(B.5) It now follows from direct substitution of the relevant terms of (B.5) into the definition of the partial correlation that: c(δ m,T ) = Cov X δ , X 2δ X 0 , X 3δ Var(X δ X 0 , X 3δ )Var(X 2δ X 0 , X 3δ ) = 1 2 cosh(θδ) . (B.6) The remaining steps are analogous to the proof of Corollary 1.
Proof of Theorem 3. The result follows easily from (7)  as m → ∞, as required. Related expressions for the relaxation time of the Ornstein-Uhlenbeck process, and for the relaxation time relating to ρ rand , follow similarly.