Library

Important concepts

Bridge.ContinuousTimeProcessType
ContinuousTimeProcess{T}

Types inheriting from the abstract type ContinuousTimeProcess{T} characterize the properties of a T-valued stochastic process, play a similar role as distribution types like Exponential in the package Distributions.

source
Bridge.SamplePathType
SamplePath{T} <: AbstractPath{T}

The struct

struct SamplePath{T}
    tt::Vector{Float64}
    yy::Vector{T}
    SamplePath{T}(tt, yy) where {T} = new(tt, yy)
end

serves as container for discretely observed ContinuousTimeProcesses and for the sample path returned by direct and approximate samplers. tt is the vector of the grid points of the observation/simulation and yy is the corresponding vector of states.

It supports getindex, setindex!, length, copy, vcat.

source
Bridge.GSamplePathType

Like VSamplePath, but with assumptions on tt and dimensionality. Planned replacement for VSamplePath

source
Base.valtypeFunction
valtype(::ContinuousTimeProcess) -> T

Returns statespace (type) of a ContinuousTimeProcess{T].

source
Bridge.outertypeFunction
outertype(P::ContinuousTimeProcess) -> T

Returns the type of outer(x), where x is a state of P

source

Ordinary differential equations and quadrature

Bridge.ODESolverType
ODESolver

Abstract (super-)type for solving methods for ordinary differential equations.

source
Bridge.solve!Function
solve!(method, Y, W, P) -> Y

Integrate with method, where P is a bridge proposal overwritingY`.

source
Bridge.R3Type
R3

Ralston (1965) update (order 3 step of the Bogacki–Shampine 1989 method) to solve $y(t + dt) - y(t) = \int_t^{t+dt} F(s, y(s)) ds$.

source
Bridge.BS3Type
BS3

Ralston (1965) update (order 3 step of the Bogacki–Shampine 1989 method) to solve $y(t + dt) - y(t) = \int_t^{t+dt} F(s, y(s)) ds$. Uses Bogacki–Shampine method to give error estimate.

source

Brownian motion

Stochastic differential equations

Bridge.aFunction
a(t, x, P::ProcessOrCoefficients)

Fallback for a(t, x, P) calling σ(t, x, P)*σ(t, x, P)'.

source
StatsBase.sampleFunction
sample(tt, P, x1=zero(T))

Sample the process P on the grid tt exactly from its transitionprob(-ability) starting in x1.

source
sample(::Thinning, T, P::InhomogPoisson) -> tt
source
StatsBase.sample!Function
sample!([::Bridge.TransitionProb], X, P, x1=zero(T))

Sample the process P on the grid X.tt exactly from its transitionprob(-ability) starting in x1 writing into X.yy.

source
Bridge.quvarFunction
quvar(X)

Computes the (realized) quadratic variation of the path X.

source
Bridge.bracketFunction
bracket(X)
bracket(X,Y)

Computes quadratic variation process of X (of X and Y).

source
Bridge.itoFunction
ito(Y, X)

Integrate a stochastic process Y with respect to a stochastic differential dX.

source
Bridge.girsanovFunction
girsanov(X::SamplePath, P::ContinuousTimeProcess, Pt::ContinuousTimeProcess)

Girsanov log likelihood $\mathrm{d}P/\mathrm{d}Pt(X)$

source
Bridge.lpFunction
lp(s, x, t, y, P)

Log-transition density, shorthand for logpdf(transitionprob(s,x,t,P),y).

source
Bridge.llikelihoodFunction
llikelihood(X::SamplePath, P::ContinuousTimeProcess)

Log-likelihood of observations X using transition density lp.

source
llikelihood(X::SamplePath, Pº::LocalGammaProcess, P::LocalGammaProcess)

Log-likelihood dPº/dP. (Up to proportionality.)

source
llikelihood(X::SamplePath, P::LocalGammaProcess)

Bridge log-likelihood with respect to reference measure P.P. (Up to proportionality.)

source
Bridge.solveFunction
solve(method::SDESolver, u, W::SamplePath, P) -> X
solve(method::SDESolver, u, W::SamplePath, (b, σ)) -> X

Solve stochastic differential equation $dX_t = b(t,X_t)dt + σ(t,X_t)dW_t$ using method in place.

Example

solve(EulerMaruyama(), 1.0, sample(0:0.1:10, Wiener()), ((t,x)->-x, (t,x)->I))
struct OU <: ContinuousTimeProcess{Float64}
    μ::Float64
end
Bridge.b(s, x, P::OU) = -P.μ*x
Bridge.σ(s, x, P::OU) = I

solve(EulerMaruyama(), 1.0, sample(0:0.1:10, Wiener()), OU(1.4))
source
solve(method::SDESolver, u, W::VSamplePath, P) -> X

Solve stochastic differential equation $dX_t = b(t,X_t)dt + σ(t,X_t)dW_t$ using method.

source
solve(method, W, P) -> Y

Integrate with method, where P is a bridge proposal from startpoint(P)toendpoint(P)`.

Examples

cs = Bridge.CSpline(tt[1], tt[end], Bridge.b(tt[1], v[1], P),  Bridge.b(tt[end], v[2], P))
P° = BridgeProp(Pσ, v), Pσ.a, cs)
W = sample(tt, Wiener())
solve(Euler(), W, P°)
source
Bridge.EulerType
EulerMaruyama() <: SDESolver

Euler-Maruyama scheme. Euler is defined as alias.

source

In place solvers

Bridge.R3!Type
R3!

Inplace Ralston (1965) update (order 3 step of the Bogacki–Shampine 1989 method) to solve $y(t + dt) - y(t) = \int_t^{t+dt} F(s, y(s)) ds$.

source
Bridge.σ!Function
σ!(t, y, Δw, tmp2, P)

Compute stochastic increment at y, $σ Δw$, modifying tmp2.

source
Bridge.b!Function
b!(t, y, tmp1, P)

Compute drift $b$ in y (without factor $Δt$, modifying tmp1.

source
Bridge.kernelr3!Function

One step for inplace Ralston (1965) update (order 3 step of the Bogacki–Shampine 1989 method) to solve $y(t + dt) - y(t) = \int_t^{t+dt} f(s, y(s)) ds$. Starting point is specified by (t,y)

f!(t,y,k) is a function that takes (t,y) and writes the result in k. ws contains 4 copies of the type of y the result is written into out which is of type y

source

Levy processes

Bridge.GammaProcessType
GammaProcess

A GammaProcess with jump rate γ and inverse jump size λ has increments Gamma(t*γ, 1/λ) and Levy measure

\[ν(x)=γ x^{-1}\exp(-λ x),\]

Here Gamma(α,θ) is the Gamma distribution in julia's parametrization with shape parameter α and scale θ.

source
Bridge.ExpCountingType
ExpCounting(λ)

Counting process with arrival times arrival(P) = Exponential(1/λ) and unit jumps.

source
Bridge.CompoundPoissonType
CompoundPoisson{T} <: LevyProcess{T}

Abstract type. For a compound Poisson process define rjumpsize(P) -> T and arrival(P) -> Distribution.

source
Bridge.nuFunction
 nu(k, P)

(Bin-wise) integral of the Levy measure $\nu(B_k)$ (sic).

source
Bridge.uniform_thinning!Function
uniform_thinning!(X, P::GammaProcess, γᵒ)

Return a Gamma process Y with new intensity γᵒ, such that X-Y has intensity γ-γᵒ and Y and X-Y are independent. In the limit $dt \to \infty$ the new Gamma process has each of is jump removed with probability γᵒ/γ. Overwrites X with Y.

source

Poisson processes

Bridge.ThinningAlgType
ThinningAlg(λmax)

Sampling method for InhomogPoisson by the 'thinning' algorithm.

Examples:

sample(ThinningAlg(λmax), T, InhomogPoisson(λ))
source
Bridge.InhomogPoissonType
InhomogPoisson(λ)

Inhomogenous Poisson process with intensity function λ(t). See also ThinningAlg.

source

Bessel processes

Bridge.BesselType
Bessel{N}(σ)

N-dimensional Bessel process with dispersion σ. Sample with

u = 0.0
t = 0:0.1:1
σ = 1.0
sample(u, t, Bridge.Bessel{3}(σ))
source
Bridge.Bessel3BridgeType
Bessel3Bridge(t, v, σ)

Bessel(3) bridge from below or above to the point v at time t, not crossing v, with dispersion σ.

source

Miscellaneous

Bridge.endpoint!Function
endpoint!(X::SamplePath, v)

Convenience functions setting the endpoint of X tov`.

source
Bridge.cumsum0Function
cumsum0(x)

Cumulative sum starting at 0 such that cumsum0(diff(x)) ≈ x.

source
Bridge.matFunction
mat(yy::Vector{SVector})

Reinterpret X or yy to an array without change in memory.

source
Bridge.CSplineType
CSpline(s, t, x, y = x, m0 = (y-x)/(t-s), m1 = (y-x)/(t-s))

Cubic spline parametrized by $f(s) = x$ and $f(t) = y$, $f'(s) = m_0$, $f'(t) = m_1$.

source
Bridge.PSDType
PSD{T}

Simple wrapper for the lower triangular Cholesky root of a positive (semi-)definite element σ.

source
Bridge.GaussianType
Gaussian(μ, Σ) -> P

Gaussian distribution with mean μand covarianceΣ. Definesrand(P)and(log-)pdf(P, x). Designed to work withNumbers,UniformScalings,StaticArraysandPSD`-matrices.

Implementation details: On Σ the functions logdet, whiten and unwhiten (or cholupper as fallback for the latter two) are called.

source
Bridge.refineFunction
refine(tt, n)

Refine range by decreasing stepsize by a factor n.

source
Bridge.quaternionFunction
quaternion(m::SMatrix{3,3})

Compute the (rotation-) quarternion of a 3x3 rotation matrix. Useful to create isodensity ellipses from spheres in GL visualizations.

source
Bridge._viridisConstant
_viridis

Color data of the Viridis map by Nathaniel J. Smith, Stefan van Der Walt, Eric Firing from https://github.com/BIDS/colormap/blob/master/colormaps.py .

source
Bridge.posteriorFunction
posterior(Val{:λ}, P::GammaProcess, U::SamplePath, prior = (0.0, 0.0))

Marginal posterior distribution of parameter λ. Interpretation of conjugate prior is "observed time, observed increment".

source

Online statistics

Online updating of the tuple state = (m, m2, n) where

m - mean(x[1:n])

m2 - sum of squares of differences from the current mean, $\textstyle\sum_{i=1}^n (x_i - \bar x_n)^2$

n - number of iterations

Bridge.mcstartFunction
mcstart(x) -> state

Create state for random chain online statitics. The entries/value of x are ignored

source
Bridge.mcnextFunction
mcnext(state, x) -> state

Update random chain online statistics when new chain value x was observed. Return new state.

source
Bridge.mcbandFunction
mcband(mc)

Compute marginal 95% coverage interval for the chain from normal approximation.

source
Bridge.mcbandmeanFunction
mcmeanband(mc)

Compute marginal confidence interval for the chain mean using normal approximation

source
Bridge.mcmarginalstatsFunction
mcmarginalstats(mcstates) -> mean, std

Compute meanand marginal standard deviationsstd` for 2d plots.

source
Bridge.OnlineStatType
stats = map(OnlineStat, (x, θ, η))

map(push!, stats, (x, θ, η))

mean.(stats)
cov.(stats)
source

Linear Processes

Bridge.PtildeType
Ptilde(cs::CSpline, σ)

Affine diffusion $dX = cs(t) dt + σdW$ with cs a cubic spline ::CSpline.

source
Bridge.LinearNoiseApprType
LinearNoiseAppr(tt, P, x, a, direction = forward)

Precursor of the linear noise approximation of P. For now no attempt is taken to add in a linearization around the deterministic path. direction can be one of :forward, :backward or :nothing. The latter corresponds to choosing β == 0.

source

Bridges

Bridge.GuidedPropType
GuidedProp

General bridge proposal process, only assuming that Pt defines H and r in the right way.

source
Bridge.GuidedBridgeType
GuidedBridge

Guided proposal process for diffusion bridge using backward recursion.

GuidedBridge(tt, P, Pt, v)

Constructor of guided proposal process for diffusion bridge of P to v on the time grid tt using guiding term derived from linear process Pt.

GuidedBridge(tt, P, Pt, V, H♢)

Guided proposal process for diffusion bridge of P to v on the time grid tt using guiding term derived from linear process Pt. Initialize using Bridge.gpupdate(H♢, V, L, Σ, v)

source
Bridge.PartialBridgeType
PartialBridge

Guided proposal process for diffusion bridge using backward recursion. PartialBridge(tt, P, Pt, L, v, Σ) Guided proposal process for a partial diffusion bridge of P to v on the time grid tt using guiding term derived from linear process Pt. Simulate with bridge!.

source
Bridge.PartialBridgeνHType
PartialBridgeνH

Guided proposal process for diffusion bridge using backward recursion.

PartialBridgeνH(tt, P, Pt, L, v,ϵ Σ)

Guided proposal process for a partial diffusion bridge of `P` to `v` on
the time grid `tt` using guiding term derived from linear process `Pt`.

PartialBridgeνH(tt, P, Pt, ν, Hend⁺)

Guided proposal process on the time grid `tt` using guiding term derived from
linear process `Pt` with backwards equation initialized at `ν, Hend⁺`.
source
Bridge.BridgePropType
BridgeProp(Target::ContinuousTimeProcess, tt, v, a, cs)

Simple bridge proposal derived from a linear process with time dependent drift given by a CSpline and constant diffusion coefficient a.

source
Bridge.MdbType
Mdb() <: SDESolver

Euler scheme with the diffusion coefficient correction of the modified diffusion bridge.

source
Bridge.VsFunction
Vs(s, T1, T2, v, P)

Time changed V for generation of U.

source
Bridge.gpV!Function

gpV!(K::SamplePath, P, KT=zero(T))

Precompute V from $(d/dt)V = BV + β$, $V_T = v$ for a guided proposal.

source
Bridge.rFunction
r(t, x, T, v, P)

Returns $r(t,x) = \operatorname{grad}_x \log p(t,x; T, v)$ where $p$ is the transition density of the process $P$.

source
Bridge.gpHinv!Function
gpHinv!(K::SamplePath, P, KT=zero(T))

Precompute $K = H^{-1}$ from $(d/dt)K = BK + KB' + a$ for a guided proposal.

source
Bridge.gpupdateFunction
gpupdate(H♢, V, L, Σ, v)
gpupdate(P, V, L, Σ, v)

Return updated H♢, V when observation v at time zero with error Σ is observed.

source

Unsorted

Bridge.compensator0Function
compensator0(kstart, P::LocalGammaProcess)

Compensator of GammaProcess approximating the LocalGammaProcess. For kstart == 1 (only choice) this is $\nu_0([b_1,\infty)$.

source
Bridge.compensatorFunction
compensator(kstart, P::LocalGammaProcess)

Compensator of LocalGammaProcess For kstart = 1, this is $\sum_{k=1}^N \nu(B_k)$, for kstart = 0, this is $\sum_{k=0}^N \nu(B_k) - C$ (where $C$ is a constant).

source
Bridge.θFunction
θ(x, P::LocalGammaProcess)

Inverse jump size compared to gamma process with same alpha and beta.

source
Bridge.softFunction
soft(t, T1, T2)

Time change mapping s in [T1, T2](U-time) tot`in[T1, T2](X`-time).

source
Bridge.tofsFunction
tofs(s, T1, T2)

Time change mapping t in [T1, T2] ($X$-time) to s in [T1, T2] (U-time).

source
Bridge.dotVsFunction
dotVs (s, T1, T2, v, P)

Time changed time derivative of V for generation of U.

source
Bridge.SDESolverType
SDESolver

Abstract (super-)type for solving methods for stochastic differential equations.

source
Bridge.IncrementsType
Increments{S<:AbstractPath{T}}

Iterator over the increments of an AbstractPath. Iterates over (i, tt[i], tt[i+1]-tt[i], yy[i+1]-y[i]).

source
Bridge.sizedtypeFunction
sizedtype(x) -> T

Return an extended type which preserves size information. Makes one(T) and zero(T) for vectors possible.

source
Bridge.piecewiseFunction
piecewise(X::SamplePath, [endtime]) -> tt, xx

If X is a jump process with piecewise constant paths and jumps in X.tt, piecewise returns coordinates path for plotting purposes. The second argument allows to choose the right endtime of the last interval.

source
Bridge.aeulerFunction
aeuler(u, s:dtmax:t, P, tau=0.5)

Adaptive Euler-Maruyama scheme from https://arxiv.org/pdf/math/0601029.pdf sampling a path from u at s to t with adaptive stepsize of 2.0^(-k)*dtmax

source
Bridge.MeanCovType
MeanCov(itr)

Iterator interface for online mean and covariance estimates. Iterates are triples mean, λ, cov/λ

Example:

c = Channel{Vector}(1)
m = Bridge.MeanCov(c)
put!(c, rand(5))
u = iterate(m)
put!(c, rand(5))
u = iterate(m, u[2])
close(c)
u[1][1]
m = Bridge.MeanCov(Channel{Vector{Float64}}(1))
u = register!(m, rand(5))
u = register!(m, rand(5), u)
close(m)
u[1][1]
source
Bridge.upsampleFunction
upsample(x, td, t)

If x is piecewise constant with jumps at td, return values of x at times t.

source
Bridge.rescaleFunction
rescale(x, a=>b, u=>v)

Linearly map the interval [a,b] to [u,v].

Example:

rescale.(x, Ref(extrema(x)))
source