# sde.sim: Simulation of stochastic differential equation In sde: Simulation and Inference for Stochastic Differential Equations

## Description

Generic interface to different methods of simulation of solutions to stochastic differential equations.

## Usage

 ```1 2 3 4 5 6 7``` ```sde.sim(t0 = 0, T = 1, X0 = 1, N = 100, delta, drift, sigma, drift.x, sigma.x, drift.xx, sigma.xx, drift.t, method = c("euler", "milstein", "KPS", "milstein2", "cdist","ozaki","shoji","EA"), alpha = 0.5, eta = 0.5, pred.corr = T, rcdist = NULL, theta = NULL, model = c("CIR", "VAS", "OU", "BS"), k1, k2, phi, max.psi = 1000, rh, A, M=1) ```

## Arguments

 `t0` time origin. `T` horizon of simulation. `X0` initial value of the process. `N` number of simulation steps. `M` number of trajectories. `delta` time step of the simulation. `drift` drift coefficient: an expression of two variables `t` and `x`. `sigma` diffusion coefficient: an expression of two variables `t` and `x`. `drift.x` partial derivative of the drift coefficient w.r.t. `x`: a function of two variables `t` and `x`. `sigma.x` partial derivative of the diffusion coefficient w.r.t. `x`: a function of two variables `t` and `x`. `drift.xx` second partial derivative of the drift coefficient w.r.t. `x`: a function of two variables `t` and `x`. `sigma.xx` second partial derivative of the diffusion coefficient w.r.t. `x`: a function of two variables `t` and `x`. `drift.t` partial derivative of the drift coefficient w.r.t. `t`: a function of two variables `t` and `x`. `method` method of simulation; see details. `alpha` weight `alpha` of the predictor-corrector scheme. `eta` weight `eta` of the predictor-corrector scheme. `pred.corr` boolean: whether to apply the predictor-correct adjustment; see details. `rcdist` a function that is a random number generator from the conditional distribution of the process; see details. `theta` vector of parameters for `cdist`; see details. `model` model from which to simulate; see details. `k1` lower bound for `psi(x)`; see details. `k2` upper bound for `psi(x)`; see details. `phi` the function `psi(x) - k1`. `max.psi` upper value of the support of `psi` to search for its maximum. `rh` the rejection function; see details. `A` `A(x)` is the integral of the `drift` between `0` and `x`.

## Details

The function returns a `ts` object of length `N+1`; i.e., `X0` and the new `N` simulated values if `M=1`. For `M>1`, an `mts` (multidimensional `ts` object) is returned, which means that `M` independent trajectories are simulated. If the initial value `X0` is not of the length `M`, the values are recycled in order to have an initial vector of the correct length. If `delta` is not specified, then `delta = (T-t0)/N`. If `delta` is specified, then `N` values of the solution of the sde are generated and the time horizon `T` is adjusted to be `N * delta`.

The function `psi` is `psi(x) = 0.5*drift(x)^2 + 0.5*drift.x(x)`.

If any of `drift.x`, `drift.xx`, `drift.t`, `sigma.x`, and `sigma.xx` are not specified, then numerical derivation is attempted when needed.

If `sigma` is not specified, it is assumed to be the constant function `1`.

The `method` of simulation can be one among: `euler`, `KPS`, `milstein`, `milstein2`, `cdist`, `EA`, `ozaki`, and `shoji`. No assumption on the coefficients or on `cdist` is checked: the user is responsible for using the right method for the process object of simulation.

The `model` is one among: `CIR`: Cox-Ingersoll-Ross, `VAS`: Vasicek, `OU` Ornstein-Uhlenbeck, `BS`: Black and Scholes. No assumption on the coefficient `theta` is checked: the user is responsible for using the right ones.

If the `method` is `cdist`, then the process is simulated according to its known conditional distribution. The random generator `rcdist` must be a function of `n`, the number of random numbers; `dt`, the time lag; `x`, the value of the process at time `t` - `dt`; and the vector of parameters `theta`.

For the exact algorithm method `EA`: if missing `k1` and `k2` as well as `A`, `rh` and `phi` are calculated numerically by the function.

## Value

 `x` returns an invisible `ts` object

## Author(s)

Stefano Maria Iacus

## References

See Chapter 2 of the text.

## Examples

 ``` 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39``` ```# Ornstein-Uhlenbeck process set.seed(123) d <- expression(-5 * x) s <- expression(3.5) sde.sim(X0=10,drift=d, sigma=s) -> X plot(X,main="Ornstein-Uhlenbeck") # Multiple trajectories of the O-U process set.seed(123) sde.sim(X0=10,drift=d, sigma=s, M=3) -> X plot(X,main="Multiple trajectories of O-U") # Cox-Ingersoll-Ross process # dXt = (6-3*Xt)*dt + 2*sqrt(Xt)*dWt set.seed(123) d <- expression( 6-3*x ) s <- expression( 2*sqrt(x) ) sde.sim(X0=10,drift=d, sigma=s) -> X plot(X,main="Cox-Ingersoll-Ross") # Cox-Ingersoll-Ross using the conditional distribution "rcCIR" set.seed(123) sde.sim(X0=10, theta=c(6, 3, 2), rcdist=rcCIR, method="cdist") -> X plot(X, main="Cox-Ingersoll-Ross") set.seed(123) sde.sim(X0=10, theta=c(6, 3, 2), model="CIR") -> X plot(X, main="Cox-Ingersoll-Ross") # Exact simulation set.seed(123) d <- expression(sin(x)) d.x <- expression(cos(x)) A <- function(x) 1-cos(x) sde.sim(method="EA", delta=1/20, X0=0, N=500, drift=d, drift.x = d.x, A=A) -> X plot(X, main="Periodic drift") ```

### Example output

```Loading required package: MASS

Attaching package: 'fda'

The following object is masked from 'package:graphics':

matplot

Attaching package: 'zoo'

The following objects are masked from 'package:base':

as.Date, as.Date.numeric

sde 2.0.15
Companion package to the book
'Simulation and Inference for Stochastic Differential Equations With R Examples'
Iacus, Springer NY, (2008)
To check the errata corrige of the book, type vignette("sde.errata")
sigma.x not provided, attempting symbolic derivation.

sigma.x not provided, attempting symbolic derivation.

sigma.x not provided, attempting symbolic derivation.

T set to = 25.000000

k1 missing, trying numerical minimization...
(k1=-0.500)

k2 missing, trying numerical maximization...
(k2=0.625)

rejection rate: 0.215
```

sde documentation built on May 31, 2017, 3:58 a.m.