# msde.fit: Estimation Of The Random Effects In Mixed Stochastic... In MsdeParEst: Parametric Estimation in Mixed-Effects Stochastic Differential Equations

## Description

Parametric estimation of the joint density of the random effects in the mixed SDE

dX_j(t)= (α_j- β_j X_j(t))dt + σ_j \ a(X_j(t)) dW_j(t),

j=1,…,M, where the (W_j(t)) are independant Wiener processes and the (X_j(t)) are observed without noise. There can be random effects either in the drift (α_j,β_j) or in the diffusion coefficient σ_j or both (α_j,β_j,σ_j).

## Usage

 1 2 3 4 msde.fit(times, X, model = c("OU", "CIR"), drift.random = c(1, 2), drift.fixed = NULL, diffusion.random = 0, diffusion.fixed = NULL, nb.mixt = 1, Niter = 10, discrete = 1, valid = 0, level = 0.05, newwindow = FALSE) 

## Arguments

 times vector of observation times X matrix of the M trajectories (each row is a trajectory with as much columns as observations) model name of the SDE: 'OU' (Ornstein-Uhlenbeck) or 'CIR' (Cox-Ingersoll-Ross) drift.random random effects in the drift: 0 if only fixed effects, 1 if one additive random effect, 2 if one multiplicative random effect or c(1,2) if 2 random effects. Default to c(1,2) drift.fixed NULL if the fixed effect(s) in the drift is (are) estimated, value of the fixed effect(s) otherwise. Default to NULL diffusion.random 1 if σ is random, 0 otherwise. Default to 0 diffusion.fixed NULL if σ is estimated (if fixed), value of σ otherwise. Default to NULL nb.mixt number of mixture components for the distribution of the random effects in the drift. Default to 1 (no mixture) Niter number of iterations for the EM algorithm if the random effects in the drift follow a mixture distribution. Default to 10 discrete 1 for using a contrast based on discrete observations, 0 otherwise. Default to 1 valid 1 if test validation, 0 otherwise. Default to 0 level alpha for the predicion intervals. Default 0.05 newwindow logical(1), if TRUE, a new window is opened for the plot. Default to FALSE

## Details

Parametric estimation of the random effects density from M independent trajectories of the SDE:

dX_j(t)= (α_j- β_j X_j(t))dt + σ_j \ a(X_j(t)) dW_j(t),

j=1,…,M, where the (W_j(t)) are independant Wiener processes and the (X_j(t)) are observed without noise.

Specification of the random effects:

The drift includes no, one or two random effects:

1. if drift.random = 0: α_j \equiv α and β_j \equiv β are fixed

2. if drift.random = 1: β_j \equiv β is fixed and α_j is random

3. if drift.random = 2: α_j \equiv α is fixed and β_j is random

4. if drift.random = c(1,2): α_j and β_j are random

The diffusion includes either a fixed effect or a random effect:

1. if diffusion.random = 0: σ_j \equiv σ is fixed

2. if diffusion.random = 1: σ_j is random

Distribution of the random effects

If there is no random effect in the diffusion (diffusion.random = 0), there is at least on random effect in the drift that follows

1. a Gaussian distribution (nb.mixt=1): α_j \sim N(μ,Ω) or β_j \sim N(μ,Ω) or (α_j,β_j) \sim N(μ,Ω),

2. or a mixture of Gaussian distributions (nb.mixt=K, K>1): α_j \sim ∑_{k=1}^{K} p_k N(μ_k,Ω_k) or β_j \sim ∑_{k=1}^{K} p_k N(μ_k,Ω_k) or (α_j,β_j) \sim ∑_{k=1}^{K} p_k N(μ_k,Ω_k), where ∑_{k=1}^{K} p_k=1.

If there is one random effect in the diffusion (diffusion.random = 1), 1/σ_j^2 \sim Γ(a,λ), and the coefficients in the drift are conditionally Gaussian: α_j|σ_j \sim N(μ,σ_j^2 Ω) or β_j|σ_j \sim N(μ,σ_j^2 Ω) or (α_j,β_j)|σ_j \sim N(μ,σ_j^2 Ω), or they are fixed α_j \equiv α, β_j \equiv β.

SDEs

Two diffusions are implemented:

1. the Ornstein-Uhlenbeck model (OU) a(X_j(t))=1

2. the Cox-Ingersoll-Ross model (CIR) a(X_j(t))=√{X_j(t)}

Estimation

• If discrete = 0, the estimation is based on the exact likelihood associated with continuous observations ([1],[3]). This is only possible if diffusion.random = 0 and σ is not estimated by maximum likelihood but empirically by means of the quadratic variations.

• If discrete = 1, the likelihood of the Euler scheme of the mixed SDE is computed and maximized for estimating all the parameters.

• If nb.mixt > 1, an EM algorithm is implemented and the number of iterations of the algorithm must be specified with Niter.

• If valid = 1, two-thirds of the sample trajectories are used for estimation, while the rest is used for validation. A plot is then provided for visual comparison between the true trajectories of the test sample and some predicted trajectories simulated under the estimated model.

## Value

 index is the vector of subscript in 1,...,M used for the estimation. Most of the time index=1:M, except for the CIR that requires positive trajectories. estimphi matrix of estimators of the drift random effects \hat{α}_j, or \hat{β}_j or (\hat{α}_j,\hat{β}_j) estimpsi2 vector of estimators of the squared diffusion random effects \hat{σ}_j^2 gridf grid of values for the plots of the random effects distribution in the drift, matrix form gridg grid of values for the plots of the random effects distribution in the diffusion, matrix form estimf estimator of the density of α_j, β_j or (α_j,β_j). Matrix form. estimg estimator of the density of σ_j^2. Matrix form. mu estimator of the mean of the random effects normal density omega estimator of the standard deviation of the random effects normal density a estimated value of the shape of the Gamma distribution lambda estimated value of the scale of the Gamma distribution sigma2 value of the diffusion coefficient if it is fixed bic BIC criterium aic AIC criterium model initial choice drift.random initial choice diffusion.random initial choice drift.fixed initial choice estim.drift.fix 1 if the fixed effects in the drift are estimated, 0 otherwise. estim.diffusion.fixed 1 if the fixed effect in the diffusion is estimated, 0 otherwise. discrete initial choice times initial choice X initial choice

For mixture distributions in the drift:

 mu estimated value of the mean at each iteration of the algorithm. Niter x N x 2 array. omega estimated value of the standard deviation at each iteration of the algorithm. Niter x N x 2 array. mixt.prop estimated value of the mixture proportions at each iteration of the algorithm. Niter x N matrix. probindi posterior component probabilites. M x N matrix.

## Author(s)

Maud Delattre and Charlotte Dion

## References

See

[1] Maximum Likelihood Estimation for Stochastic Differential Equations with Random Effects, Delattre, M., Genon-Catalot, V. and Samson, A. Scandinavian Journal of Statistics 40(2) 2012 322-343

[2] Estimation of population parameters in stochastic differential equations with random effects in the diffusion coefficient, Delattre, M., Genon-Catalot, V. and Samson, A. ESAIM:PS 19 2015 671-688

[3] Mixtures of stochastic differential equations with random effects: application to data clustering, Delattre, M., Genon-Catalot, V. and Samson, A. Journal of Statistical Planning and Inference 173 2016 109-124

[4] Parametric inference for discrete observations of diffusion processes with mixed effects, Delattre, M., Genon-Catalot, V. and Laredo, C. hal-01332630 2016

[5] Estimation of the joint distribution of random effects for a discretely observed diffusion with random effects, Delattre, M., Genon-Catalot, V. and Laredo, C. hal-01446063 2017

## 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 40 41 42 43 44 45 46 47 48 49 50 51 # Example 1 : One random effect in the drift and one random effect in the diffusion sim <- msde.sim(M = 25, T = 1, N = 1000, model = 'OU', drift.random = 2, drift.param = c(0,0.5,0.5), diffusion.random = 1, diffusion.param = c(8,1/2)) res <- msde.fit(times = sim$times, X = sim$X, model = 'OU', drift.random = 2, diffusion.random = 1) summary(res) plot(res) ## Not run: # Example 2 : one mixture of two random effects in the drift, and one fixed effect in # the diffusion coefficient sim <- msde.sim(M = 100, T = 5, N = 5000, model = 'OU', drift.random = c(1,2), diffusion.random = 0, drift.param = matrix(c(0.5,1.8,0.25,0.25,1,2,0.25,0.25),nrow=2,byrow=FALSE), diffusion.param = 0.1, nb.mixt = 2, mixt.prop = c(0.5,0.5)) # -- Estimation without validation res <- msde.fit(times = sim$times, X = sim$X, model = 'OU', drift.random = c(1,2), nb.mixt=2, Niter = 10) summary(res) plot(res) # -- Estimation with prediction res.valid <- msde.fit(times = sim$times, X = sim$X, model = 'OU', drift.random = c(1,2), nb.mixt=2, Niter = 10, valid = 1) summary(res.valid) plot(res.valid) # Example 3 : CIR with one random effect in the drift and one random effect in the diffusion # coefficient sim <- msde.sim(M = 100, T = 5, N = 5000, model = 'CIR', drift.random = 2, diffusion.random = 1, drift.param = c(4,1,0.1), diffusion.param = c(8,0.5), X0 = 1) res <- msde.fit(times = sim$times, X = sim$X, model = 'CIR', drift.random = 2, diffusion.random = 1) summary(res) # Further examples can be found in the section dedicated to neuronal.data. ## End(Not run) 

MsdeParEst documentation built on May 2, 2019, 3:01 a.m.