rTrajLangevin: Simulation of trajectories of a Langevin diffusion

In egarpor/sdetorus: Statistical Tools for Toroidal Diffusions

Simulation of trajectories of a Langevin diffusion

Description

Simulation of an arbitrary Langevin diffusion in dimension p by subsampling a fine trajectory obtained by the Euler discretization.

Usage

rTrajLangevin(x0, drift, SigDif, N = 100, delta = 0.01, NFine = ceiling(N
  * delta/deltaFine), deltaFine = min(delta/100, 0.001), circular = TRUE,
  ...)

Arguments

x0	vector of length p giving the initial point.
drift	drift for the diffusion.
SigDif	matrix of size c(p, p) giving the infinitesimal (constant) covariance matrix of the diffusion.
N	number of discretization steps in the resulting trajectory.
delta	discretization step.
NFine	number of discretization steps for the fine trajectory. Must be larger than N.
deltaFine	discretization step for the fine trajectory. Must be smaller than delta.
circular	whether to wrap the resulting trajectory to [-\pi,\pi)^p.
...	parameters to be passed to drift.

Details

The fine trajectory is subsampled using the indexes seq(1, NFine + 1, by = NFine / N).

Value

A vector of length N + 1 containing x0 in the first entry and the discretized trajectory.

Examples

isRStudio <- identical(.Platform$GUI, "RStudio")
if (isRStudio) {
  # 1D
  manipulate::manipulate({
    x <- seq(0, N * delta, by = delta)
    plot(x, x, ylim = c(-pi, pi), type = "n",
         ylab = expression(X[t]), xlab = "t")
    linesCirc(x, rTrajLangevin(x0 = 0, drift = driftJp, SigDif = sigma,
                               alpha = alpha, mu = 0, psi = psi, N = N,
                               delta = 0.01))
    }, delta = manipulate::slider(0.01, 5.01, step = 0.1),
    N = manipulate::slider(10, 500, step = 10, initial = 200),
    alpha = manipulate::slider(0.01, 5, step = 0.1, initial = 1),
    psi = manipulate::slider(-2, 2, step = 0.1, initial = 1),
    sigma = manipulate::slider(0.01, 5, step = 0.1, initial = 1))

  # 2D
  samp <- rTrajLangevin(x0 = c(0, 0), drift = driftMvm, alpha = c(1, 1),
                        mu = c(2, -1), A = diag(rep(0, 2)),
                        SigDif = diag(rep(1, 2)), N = 1000, delta = 0.1)
  plot(samp, xlim = c(-pi, pi), ylim = c(-pi, pi), pch = 19, cex = 0.25,
       xlab = expression(X[t]), ylab = expression(Y[t]), col = rainbow(1000))
  linesTorus(samp[, 1], samp[, 2], col = rainbow(1000))
}

egarpor/sdetorus documentation built on March 4, 2024, 1:23 a.m.