LinProc

Introduction

This module covers

  • The simulation of multivariate diffusion processes with a simple Euler scheme
  • the simulation of (Vector) Ornstein–Uhlenbeck processes
  • the simulation of Ornstein–Uhlenbeck bridges
  • mean and covariance functions and transition density of Ornstein–Uhlenbeck processes
  • the Monte Carlo simulation of Diffusion bridges

A multivariate diffusion process is the solution to the stochastic differential equation (SDE)

\(dX_t = b(t, X_t) \sigma(t, X_t) d W_t\)

where W_t is a d-dimension Wiener process, b is a vector valued drift function and a = sigma sigma' the diffusion matrix. The function

eulerv(t0, u, b(s,x), sigma(s,x), Dt, DW::Matrix)

implements the multivariate euler scheme, starting in the point u the same conventions as in module Diffusion apply.

A vector linear process (Ornstein–Uhlenbeck process) is of the special form

\(dX_t = B X_t + \beta + \sigma d W_t\)

where B is a stable matrix and beta a vector, A = sigma sigma'

A conditional vector linear processes (Ornstein–Uhlenbeck bridges so to say) ending at time T in point v are given by

\(dX^\star_t = B X_t + \beta + A r(s, X^\star_t) + \sigma d W_t\)

where \(r(t,x) = \operatorname{grad}_x \log p(t,x; T, v)\) and p is the transition density of the corresponding linear process X.

The parameter lambda is the solution to the Lyapunov equation B lambda + lambda B' = -A, see module Lyap,

lambda = lyap(B', -A)

If B = 0, lambda = lyap(B',-A) is not defined, provide lambda = inv(a) as argument to the functions instead.

Reference

LinProc.mu(h, x, B, beta)

Expectation \(E_x(X_{t})\)

LinProc.K(h, B, lambda)

Covariance matrix \(Cov(X_{t}, X_{t + h})\)

LinProc.r(h, x, v, B, beta, lambda)

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

LinProc.H(h, B, lambda)

Negative Hessian of \(\log p(t,x; T, v)\) as a function of x.

LinProc.bstar(T, v, b, beta, a, lambda)

Returns the drift function of a vector linear process bridge which end at time T in point v.

LinProc.bcirc(T, v, b, beta, a, lambda)

Drift for guided proposal derived from a vector linear process bridge which end at time T in point v.

LinProc.llikeliXcirc(t, T, Xcirc, b, a, B, beta, lambda)

Loglikelihood (log weights) of Xcirc with respect to Xstar.

t, T – timespan Xcirc – bridge proposal (drift Bcirc and diffusion coefficient sigma) b, sigma – diffusion coefficient sigma target B, beta – drift b(x) = Bx + beta of Xtilde lambda – solution of the lyapunov equation for Xtilde
LinProc.lp(h, x, y, b, beta, lambda)

Returns \(log p(t,x; T, y)\), the log transition density of the linear process, h = T - t

LinProc.sample_p(h, x, b, beta, lambda)

Samples from the transition density of the linear process, h = T - t.

LinProc.linexact(u, B, beta, lambda, dt)

Simulate linear process starting in u on a discrete grid dt from its transition probability, corresponding to drift parameters B, beta and Lyapunov matrix lambda.

LinProc.linll(X, B, beta, lambda, dt)

Compute log likelihood evaluated in B, beta and Lyapunov matrix lambda for a observed linear process on a discrete grid dt from its transition density.

LinProc.lp0(h, x, y, mu, gamma)

Returns \(log p(t,x; T, y)\), the log transition density of a Brownian motion with drift mu and diffusion a=inv(gamma), h = T - t

LinProc.sample_p0(h, x, mu, l)

Samples from the transition density a affine Brownian motion. Takes the Cholesky factor as argument.

l = chol(a)
LinProc.eulerv(t0, u, v, b(s, x), sigma(s, x), Dt, DW::Matrix)
LinProc.eulerv(t0, u, b, sigma, dt, dw::Matrix) = eulerv(t0, u, NaN, b, sigma, dt, dw::Matrix)

Multivariate euler scheme, starting in u, fixing X[N] = v if v!=NaN (this makes sense, if b pulls X towards v.). dw – Wiener differential with n values in the the interval [t0,sum(dt)] sampled at timepoints t0+Dt[1], t0 + Dt[1] + Dt[2], ... b, sigma – drift and diffusion coefficient.

Example:

Dt = diff(linspace(0., T, N)) DW = randn(2, N-1) .* sqrt(dt) dt = Dt[1] yy = euler(0.0, u, b, sigma, Dt, DW)
LinProc.stable(Y, d, ep)

Return real stable d-dim matrix with real eigenvalues smaller than ep parametrized with a vector of length d*d,

For maximum likelihood estimation we need to search the maximum over all stable matrices. These are matrices with eigenvalues with strictly negative real parts. We obtain a dxd stable matrix as difference of a antisymmetric matrix and a positive definite matrix.