mixedsde.sim: Simulation Of A Mixed Stochastic Differential Equation In mixedsde: Estimation Methods for Stochastic Differential Mixed Effects Models

Description

Simulation of M independent trajectories of a mixed stochastic differential equation (SDE) with linear drift and two random effects (α_j, β_j) dX_j(t)= (α_j- β_j X_(t))dt + σ a(X_j(t)) dW_j(t), for j=1, ..., M.

Usage

 1 2 3 mixedsde.sim(M, T, N = 100, model, random, fixed = 0, density.phi, param, sigma, t0 = 0, X0 = 0.01, invariant = 0, delta = T/N, op.plot = 0, add.plot = FALSE)

Arguments

 M number of trajectories T horizon of simulation. N number of simulation steps, default Tx100. model name of the SDE: 'OU' (Ornstein-Uhlenbeck) or 'CIR' (Cox-Ingersoll-Ross). random random effects in the drift: 1 if one additive random effect, 2 if one multiplicative random effect or c(1,2) if 2 random effects. fixed fixed effects in the drift: value of the fixed effect when there is only one random effect, 0 otherwise. If random =2, fixed can be 0 but β has to be a non negative random variable for the estimation. density.phi name of the density of the random effects. param vector of parameters of the distribution of the two random effects. sigma diffusion parameter t0 time origin, default 0. X0 initial value of the process, default X0=0. invariant 1 if the initial value is simulated from the invariant distribution, default 0.01 and X0 is fixed. delta time step of the simulation (T/N). op.plot 1 if a plot of the trajectories is required, default 0. add.plot 1 for add trajectories to an existing plot

Details

Simulation of M independent trajectories of the SDE (the Brownian motions Wj are independent), with linear drift. Two diffusions are implemented, with one or two random effects:

Ornstein-Uhlenbeck model (OU)

If random = 1, β is a fixed effect: dX_j(t)= (α_j- β X_j(t))dt + σ dW_j(t)

If random = 2, α is a fixed effect: dX_j(t)= (α - β_j X_j(t))dt + σ dW_j(t)

If random = c(1,2), dX_j(t)= (α_j- β_j X_j(t))dt + σ dW_j(t)

Cox-Ingersoll-Ross model (CIR)

If random = 1, β is a fixed effect: dX_j(t)= (α_j- β X_j(t))dt + σ √{X_j(t)} dW_j(t)

If random = 2, α is a fixed effect: dX_j(t)= (α - β_j X_j(t))dt + σ √{X_j(t)} dW_j(t)

If random = c(1,2), dX_j(t)= (α_j- β_j X_j(t))dt + σ √{X_j(t)} dW_j(t)

The initial value of each trajectory can be simulated from the invariant distribution of the process: Normal distribution with mean α/β and variance σ^2/(2 β) for the OU, a gamma distribution Γ(2α/σ^2, σ^2/(2β)) for the C-I-R model.

Density of the random effects

Several densities are implemented for the random effects, depending on the number of random effects.

If two random effects, choice between

'normalnormal': Normal distributions for both α β and param=c(mean_α, sd_α, mean_β, sd_β)

'gammagamma': Gamma distributions for both α β and param=c(shape_α, scale_α, shape_β, scale_β)

'gammainvgamma': Gamma for α, Inverse Gamma for β and param=c(shape_α, scale_α, shape_β, scale_β)

'normalgamma': Normal for α, Gamma for β and param=c(mean_α, sd_α, shape_β, scale_β)

'normalinvgamma': Normal for α, Inverse Gamma for β and param=c(mean_α, sd_α, shape_β, scale_β)

'gammagamma2': Gamma +2 * σ^2 for α, Gamma + 1 for β and param=c(shape_α, scale_α, shape_β, scale_β)

'gammainvgamma2': Gamma +2 * σ^2 for α, Inverse Gamma for β and param=c(shape_α, scale_α, shape_β, scale_β)

If only α is random, choice between

'normal': Normal distribution with param=c(mean, sd)

lognormal': logNormal distribution with param=c(mean, sd)

'mixture.normal': mixture of normal distributions p N(μ1,σ1^2) + (1-p)N(μ2, σ2^2) with param=c(p, μ1, σ1, μ2, σ2)

'gamma': Gamma distribution with param=c(shape, scale)

'mixture.gamma': mixture of Gamma distribution p Γ(shape1,scale1) + (1-p)Γ(shape2,scale2) with param=c(p, shape1, scale1, shape2, scale2)

'gamma2': Gamma distribution +2 * σ^2 with param=c(shape, scale)

'mixed.gamma2': mixture of Gamma distribution p Γ(shape1,scale1) + (1-p) Γ(shape2,scale2) + +2 * σ^2 with param=c(p, shape1, scale1, shape2, scale2)

If only β is random, choice between 'normal': Normal distribution with param=c(mean, sd)

'gamma': Gamma distribution with param=c(shape, scale)

'mixture.gamma': mixture of Gamma distribution p Γ(shape1,scale1) + (1-p) Γ(shape2,scale2) with param=c(p, shape1, scale1, shape2, scale2)

Value

 X matrix (M x (N+1)) of the M trajectories. phi vector (or matrix) of the M simulated random effects.

References

This function mixedsde.sim is based on the package sde, function sde.sim. See Simulation and Inference for stochastic differential equation, S.Iacus, Springer Series in Statistics 2008 Chapter 2

Examples

 1 2 3 4 5 6 7 #Simulation of 5 trajectories of the OU SDE with random =1, and a Gamma distribution. simuOU <- mixedsde.sim(M=5, T=10,N=1000,model='OU', random=1,fixed=0.5, density.phi='gamma', param=c(1.8, 0.8) , sigma=0.1,op.plot=1) X <- simuOU\$X ; phi <- simuOU\$phi hist(phi)

Example output  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")
be careful, X0 and invariant are missing thus the initial value X0=0.01 is used

T set to = 10.000000

T set to = 10.000000

T set to = 10.000000

T set to = 10.000000

T set to = 10.000000

mixedsde documentation built on May 1, 2019, 7:33 p.m.