[Yuima-commits] r243 - pkg/yuimadocs/inst/doc/JSS
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Thu Jul 18 12:42:26 CEST 2013
Author: iacus
Date: 2013-07-18 12:42:25 +0200 (Thu, 18 Jul 2013)
New Revision: 243
Modified:
pkg/yuimadocs/inst/doc/JSS/article-new.Rnw
pkg/yuimadocs/inst/doc/JSS/bibliography.bib
Log:
final version of submission
Modified: pkg/yuimadocs/inst/doc/JSS/article-new.Rnw
===================================================================
--- pkg/yuimadocs/inst/doc/JSS/article-new.Rnw 2013-04-18 01:39:08 UTC (rev 242)
+++ pkg/yuimadocs/inst/doc/JSS/article-new.Rnw 2013-07-18 10:42:25 UTC (rev 243)
@@ -70,7 +70,7 @@
%% almost as usual
\author{Alexandre Brouste\\University of Le Mans \And
Masaaki Fukasawa\\Osaka University\And
- Hideitsu Hino\\Waseda University\And
+ Hideitsu Hino\\University of Tsukuba\And
Stefano M. Iacus\\University of Milan \AND
Kengo Kamatani\\University of Tokyo \And
Yuta Koike\\University of Tokyo\And
@@ -169,12 +169,12 @@
\section{Introduction}
-The plan of the YUIMA Project is to construct the bases for a complete environment for simulation and inference for stochastic processes via an \proglang{R} \citep{ERRE} package called \pkg{yuima}.
+The plan of the YUIMA\footnote{YUIMA is both the acronym for \textit{Yoshida-Uchida-Iacus-Masuda-Andothers} but also the name of an important character in Buddhism religion (\url{http://www.kyohaku.go.jp/eng/syuzou/meihin/kaiga/chuugoku/item01.html}) whose approach to problems fits well this project.} Project is to construct the bases for a complete environment for simulation and inference for stochastic processes via an \proglang{R} \citep{ERRE} package called \pkg{yuima}.
% The package \pkg{yuima} provides an object-oriented programming environment for simulation and statistical inference for stochastic processes by \proglang{R}.
The \pkg{yuima} package adopts the S4 system of classes and methods \citep{chambers98}.
Under this framework,
-the \pkg{yuima} package also supplies various functions
+the \pkg{yuima} package also supplies various functions
to carry out simulation and statistical analysis.
Both categories of procedures may depend on each other.
Statistical inference often requires a simulation technique
@@ -216,7 +216,7 @@
The \pkg{yuima} package is intended to offer the basic infrastructure on which complex models and inference procedures can be built on. This paper explains the design of the \pkg{yuima} package and illustrates some examples of applications.
-The paper is organized as follows. Section \ref{sec2} is an overview about the package. Section \ref{sec3} describes the way in which models are specified in \pkg{yuima}. Section \ref{sec4} explains asymptotic expansion methods. Section \ref{sec5} is a review of basic inference procedures.
+The paper is organized as follows. Section \ref{sec2} is an overview about the package. Section \ref{sec3} describes the way in which models are specified in \pkg{yuima}. Section \ref{sec:simul} describes how to simulate \code{yuima} models. Section \ref{sec4} explains asymptotic expansion methods. Section \ref{sec5} is a review of basic inference procedures.
Finally, Section \ref{sec6} gives additional details and the roadmap of the YUIMA Project.
@@ -253,23 +253,23 @@
\label{fig:classes}
\end{figure}
The different slots do not need to be all present at the same time. For example, in case one wants to simulate a stochastic process, only the slots \code{model} and \code{sampling} should be present, while the slot \code{data} will be filled by the simulator.
-We discuss in details the different objects separately in the next sections.
+We discuss in detail the different objects separately in the next sections.
-The general idea of the \pkg{yuima} package is to separate into different subclass objects the statistical model, the data and the statistical methods. As will be explained with several examples, the user may give a mathematical description of the statistical model with \code{setModel} which prepares a \code{yuima.model} object by filling the appropriate slots. If the aim is the simulation of the solution of the stochastic differential equation specified in the \code{yuima.model} object then, using the method \code{simulate}, it is possible to obtain one trajectory of the process. As an output, a \code{yuima} object is created which contains the original model specified in the \code{yuima.model} object in the slot named \code{model} and two additional slots named \code{data}, for the simulated data, and \code{sampling} which contains the description of the simulation scheme used as well as other informations. The details of \code{simulate} will be explained in Section \ref{sec:simul} along with the use of method \code{setSampling} which allows to specify the type of sampling scheme to be used by the \code{simulate} method.
+The general idea of the \pkg{yuima} package is to separate into different subclass objects the statistical model, the data and the statistical methods. As will be explained with several examples, the user may give a mathematical description of the statistical model with \code{setModel} which prepares a \code{yuima.model} object by filling the appropriate slots. If the aim is the simulation of the solution of the stochastic differential equation specified in the \code{yuima.model} object then, using the method \code{simulate}, it is possible to obtain one trajectory of the process. As an output, a \code{yuima} object (by ``\code{yuima} object'' we mean a, possibly incomplete, object of class \pkg{yuima}) is created which contains the original model specified in the \code{yuima.model} object in the slot named \code{model} and two additional slots named \code{data}, for the simulated data, and \code{sampling} which contains the description of the simulation scheme used as well as other informations. The details of \code{simulate} will be explained in Section \ref{sec:simul} along with the use of method \code{setSampling} which allows to specify the type of sampling scheme to be used by the \code{simulate} method.
But \code{yuima} object may contain the slot \code{data} non only as the outcome of \code{simulate} but also because the user decides to analyse its own data. In this case the method \code{setData} is used to transform most types of \proglang{R} time series object into a a proper \code{yuima.data} object.
-When the slots \code{data} and \code{model} are available, many other methods can be used to perform statistical analysis on these sde's models. These methods will be discussed in Section \ref{sec5}.
-Further, functionals of stochastic differential equations can be defined using the \code{setFunctional} method and evaluated using asymptotc expansion methods as explained in Section \ref{sec4}. The \code{setFunctional} method creates a \code{yuima.functional} object which is included along with a \code{yuima.model} into a \code{yuima} object in order to be used for the evaluation of its asymptotic expansion.
+When the slots \code{data} and \code{model} are available, many other methods can be used to perform statistical analysis on these SDE's (stochastic differential equation) models. These methods will be discussed in Section \ref{sec5}.
+Further, functionals of stochastic differential equations can be defined using the \code{setFunctional} method and evaluated using asymptotic expansion methods as explained in Section \ref{sec4}. The \code{setFunctional} method creates a \code{yuima.functional} object which is included along with a \code{yuima.model} into a \code{yuima} object in order to be used for the evaluation of its asymptotic expansion.
\subsection{The \code{yuima.model} Class}\label{sec:model}
At present, in \pkg{yuima} three main classes of stochastic differential equations can be easily specified. Here we present a brief
-overview of these models as they will be described in details in Section \ref{sec3}, but this allow to introduce an overall view of the slots of the \code{yuima.model} class.
+overview of these models as they will be described in detail in Section \ref{sec3}, but this allows to introduce an overall view of the slots of the \code{yuima.model} class.
In \pkg{yuima} one can describe three main families of stochastic processes at present. These models can be one or multidimensional and eventually described as parametric models. Let $X_0=x_0$ be the initial value of the process, then, the main classes are as follows:
\begin{itemize}
\item diffusion models described as
$$ \de X_t=a(t,X_t,\theta)dt + b(t,X_t,\theta)\de W_t, \quad X_0=x_0$$
where $W_t$ is a standard Brownian motion;
-\item sde's driven by fractional Gaussian noise, with $H$ the Hurst parameter, described as
+\item SDE's driven by fractional Gaussian noise, with $H$ the Hurst parameter, described as
$$ \de X_t=a(t,X_t,\theta)dt + b(t,X_t,\theta)\de W_t^{H}, \quad X_0=x_0;$$
\item diffusion process with jumps and L\'evy processes solution to
$$
@@ -279,8 +279,8 @@
\end{aligned}, \quad X_0=x_0;
$$
\end{itemize}
- The functions $a(\cdot)$, $b(\cdot)$ and $c(\cdot)$ may have a different number of arguments. For example, if the model is homogeneous in time and drift and diffusion coefficients are entirely specified, then we will use the notaion $a(x)$ and $b(x)$ and describe the diffusion model simply as $\displaystyle \de X_t=a(X_t)dt + b(X_t)\de W_t$. And so forth. Detailed hypotheses and regularity conditions on the coefficients $a(\cdot)$, $b(\cdot)$ and $c(\cdot)$ for each class of models will be given in the next sections. Nevertheless, it is important to remark that these notations only matter the mathematical description of the model while each coefficient is passed to \pkg{yuima} methods as \proglang{R} mathematical expressions. It means that, for example, $a(t, X_t, \theta) = t\cdot \sqrt{\theta X_t}$ will be passed as \code{t*sqrt(x*theta)}, therefore, the order of the arguments is not relevant to \proglang{R} as well as its mathematical description as long as it is consistent through each specific section. Furhter, \pkg{yuima} is able to accept any user-specified notation for the state variable $x$ (for $X_t$) and the time variable $t$ and the remaining term of an \proglang{R} expression will be interpreted as parameter as will explained in Section \ref{sec:diff}.
-We are now able to give an overview of the mian slots of the most important class of the \pkg{yuima} package.
+ The functions $a(\cdot)$, $b(\cdot)$ and $c(\cdot)$ may have a different number of arguments. For example, if the model is homogeneous in time and drift and diffusion coefficients are entirely specified, then we will use the notaion $a(x)$ and $b(x)$ and describe the diffusion model simply as $\displaystyle \de X_t=a(X_t)dt + b(X_t)\de W_t$. And so forth. Detailed hypotheses and regularity conditions on the coefficients $a(\cdot)$, $b(\cdot)$ and $c(\cdot)$ for each class of models will be given in the next sections. Nevertheless, it is important to remark that these notations only matter to the mathematical description of the model while each coefficient is passed to \pkg{yuima} methods as \proglang{R} mathematical expressions. This means that, for example, $a(t, X_t, \theta) = t\cdot \sqrt{\theta X_t}$ will be passed as \code{t*sqrt(x*theta)}, therefore, the order of the arguments is not relevant to \proglang{R} as well as its mathematical description as long as it is consistent through each specific section. Furhter, \pkg{yuima} is able to accept any user-specified notation for the state variable $x$ (for $X_t$) and the time variable $t$ and the remaining term of an \proglang{R} expression will be interpreted as parameter as will be explained in Section \ref{sec:diff}.
+We are now able to give an overview of the main slots of the most important class of the \pkg{yuima} package.
The \code{yuima.model} class contains information about the stochastic differential equation of interest. The constructor function \code{setModel} is
used to give a mathematical description of the stochastic differential equation.
@@ -298,7 +298,7 @@
coefficient matrix relative to first noise;
\item \code{hurst} is the Hurst index of the fractional Brownian motion, by default \code{0.5} meaning a standard Brownian motion. More details will be given in Section \ref{sec:fgn};
\item \code{parameter}, which is a short name for ``parameters'',
- is a list with the following entries (explained details in Section \ref{sec:par}):
+ is a list with the following entries (more details in Section \ref{sec:par}):
\begin{itemize}
\item \code{all} contains the names of all the parameters found
in the diffusion and drift coefficient;
@@ -311,26 +311,26 @@
\item \code{measure} specifies the measure of the L\'evy noise (see Section \ref{sec:jump});
\item \code{measure.type} a switch to specify how the L\'evy measure is described (see Section \ref{sec:jump});
\item \code{state.variable} and \code{time.variable}, by default,
- are assumed to be \code{x} and \code{t} but the user can freely choose them and they matter the right hand-side of the equation of the sde.
+ are assumed to be \code{x} and \code{t} but the user can freely choose them and they matter to the right hand-side of the equation of the SDE.
The \code{yuima.model} class assumes that the user either uses default
- names for \code{state.variable} and \code{time.variable} variables or specify his own names. All
- the rest of the symbols are considered parameters and distributed accordingly
+ names for \code{state.variable} and \code{time.variable} variables or specifies his own names. All
+ other symbols are considered parameters and distributed accordingly
in the \code{parameter} slot. Example of use will be given in Section \ref{sec:diff};
-\item \code{jump.variable} the name of the variable used in the description of the L\'evy component (see Section \ref{sec:jump});
+\item \code{jump.variable} indicates the name of the variable used in the description of the L\'evy component (see Section \ref{sec:jump});
\item \code{solve.variable} contains a vector of variable names, each element corresponds to the
name of the solution variable (left-hand-side) of each equation in the model, in the corresponding order. An example of use can be found in Section \ref{sec:multi}.
\item \code{noise.number} indicates the number of sources of noise.
-\item \code{xinit} initial value of the stochastic differential equation;
+\item \code{xinit} is the initial value of the stochastic differential equation;
\item \code{equation.number} represents the number of equations, i.e., the number of one dimensional
stochastic differential equations.
\item \code{dimension} reports the dimensions of the parameter space. It is a list of the
same length of \code{parameter} with the same names.
-\item \code{J.flag} for internal use only.
+\item \code{J.flag} is for internal use only.
\end{itemize}
-As seen in the above, the parameter space is accurately described internally in a \code{yuima} object because in inference for stochastic differential equations, estimators of different parameters have different properties. Usually, the rate of convergence for estimators in the diffusion coefficient are similar to the ones in the i.i.d. sampling while estimators of parameters in the drift coefficient are slower or, in some cases, not even consistent. The \pkg{yuima} always tries to implement the best statistical inference for the given model under the observed sampling scheme.
+As seen in the above, the parameter space is accurately described internally in a \code{yuima} object because in inference for stochastic differential equations, estimators of different parameters have different properties. Usually, the rate of convergence for estimators in the diffusion coefficient are similar to the ones in the i.i.d. (independent and idensitcally distributed) sampling while estimators of parameters in the drift coefficient are slower or, in some cases, not even consistent. The \pkg{yuima} always tries to implement the best statistical inference for the given model under the observed sampling scheme.
\section{Model Specification}\label{sec3}
-In order to show how general the approach is in the \pkg{yuima} package we present some examples. Throughout this section we assume that all the stochastic differential equations exists while in Section \ref{sec5} we will give precise regularity conditions needed to have a properly defined statistical model.
+In order to show how general the approach of the \pkg{yuima} package is, we present some examples. Throughout this section we assume that all the stochastic differential equations exist while in Section \ref{sec5} we will give regularity conditions needed to have a properly defined statistical model.
\subsection{One Dimensional Diffusion Processes}\label{sec:diff}
Assume that we want to describe the following stochastic differential equation
@@ -344,15 +344,15 @@
<<echo=TRUE, print=FALSE,results=hide>>=
mod1 <- setModel(drift = "-3*x", diffusion = "1/(1+x^2)")
@
-By default, \pkg{yuima} assumes that the state variable is \code{x} and the time variable is \code{t} and the solve variable is again \code{x}. Notice that the left hand-side of the equation is implicit, this is why \code{yuima.model} has the slot \code{solve.variable.
-The user should not we worried about the warning raised by \pkg{yuima} at this stage, as this is just to inform the user about the implicit assumption on the solution variable.}
+By default, \pkg{yuima} assumes that the state variable is \code{x} and the time variable is \code{t} and the solve variable is again \code{x}. Notice that the left hand-side of the equation is implicit, this is why \code{yuima.model} has the slot \code{solve.variable}.
+The user should not be worried about the warning raised by \pkg{yuima} at this stage, as this is just to inform her on the implicit assumption on the solution variable.
At this point, the package fills the proper slots of the \code{yuima} object
<<>>=
str(mod1)
@
From the above, it is possible to see that the jump coefficient is void and the Hurst parameter is set to 0.5, because this is a model where the driving process is the standard Brownian motion, i.e. a fractional Brownian motion if Hurst index $H=\frac12$.
-Now, with \code{mod1} in hands, it is very easy to simulate a trajectory by Euler-Maruyama scheme of the process as follows
+Now, with \code{mod1} in hands, it is extremely easy to simulate a trajectory by Euler-Maruyama scheme of the process as follows
<<echo=TRUE, print=FALSE,fig=TRUE,width=9,height=4,results=hide>>=
set.seed(123)
X <- simulate(mod1)
@@ -368,7 +368,7 @@
\subsection{User Specified State and Time Variables}
-Suppose now that the user wants to specify her own model using a prescribed notation, e.g., some sdes like
+Suppose now that the user wants to specify her own model using a prescribed notation, e.g., some SDE's like
$$\de Y_s = -3 s Y_s \de s + \frac{1}{1+Y_s^2}\de W_s,$$
where $a(s,y) = -3sy$ and $b(y) = 1/(1+y^2)$.
Then this model can be described in \pkg{yuima} as follows
@@ -416,7 +416,7 @@
\left(\begin{array}{ccc}1 & 0 & X_{2,t} \\X_{1,t} & 3 & 0\end{array}\right)
\left(\begin{array}{c}\de W_{1,t} \\ \de W_{2,t} \\ \de W_{3,t}\end{array}\right)
$$
-For this system it now necessary to instruct \pkg{yuima} about the state variable on bot the left-hand side of the equation and the right-hand side of the equation, i.e. the \code{solve.variable}.
+For this system it is now necessary to instruct \pkg{yuima} about the state variable on both the left-hand side of the equation and the right-hand side of the equation, i.e. there is the need to specify also the \code{solve.variable} for the left hand side of the SDE:
<<echo=TRUE, print=FALSE,results=hide>>=
sol <- c("x1","x2") # variable for numerical solution
a <- c("-3*x1","-x1-2*x2") # drift vector
@@ -431,25 +431,35 @@
plot(X, plot.type="single",lty=1:2)
@
+Notice that the role of \code{solve.variable} is essential because it allows to specify the left hand side of an SDE equation. For example, if one wants to specify this model instead of the previous one
+
+$$
+\begin{aligned}
+\de X_{2,t} &= -3 X_{1,t} \de t + \de W_{1,t} + X_{2,t} \de W_{3,t}\\
+\de X_{1,t} &= -(X_{1,t} + 2 X_{2,t}) \de t + X_{1,t} \de W_{1,t} + 3 \de W_{2,t}
+\end{aligned}
+$$
+the \code{solve.variable} argument should be specified as \code{solve.variable=c("x2","x1")} in place of \code{solve.variable=c("x1","x2")}, all the rest being the same as in model \code{mod3}.
+
\noindent
But it is also possible to specify more complex models like the following
- \begin{equation} \label{sabr}
+$$
\begin{cases}
& \mbox{d} X_{1,t} = X_{2,t} \left| X_{1,t} \right|^{2/3} \mbox{d}W_{1,t}, \\
& \mbox{d}X_{2,t} = g(t)\mbox{d}X_{3,t}, \\
& \mbox{d}X_{3,t} = X_{3,t}( \mu \mbox{d}t + \sigma (\rho \mbox{d}W_{1,t} + \sqrt{1-\rho^2}
\mbox{d}W_{2,t}))
\end{cases}
- \end{equation}
+$$
$$ (X_{1,0}, X_{2,0}, X_{3,0}) = (1.0,0.1,1.0)$$
with $\mu = 0.1, \sigma = 0.2, \rho = -0.7 $ and $g(t) = 0.4 + (0.1 + 0.2t) e^{-2t}$, for
example, where $W = (W_1, W_2)$ is a 2-dimensional standard Brownian motion. In order to pass this model to \pkg{yuima} we need to rewrite it in matrix form.
The solution $X = (X_1, X_2, X_3)$ takes values on $\mathbb{R}_+^3$.
This is a stochastic integral equation defined as
- \begin{equation}
+ $$
X_t = X_0 + \int_0^t a(s,X_s) \mbox{d}s + \int_0^t b(s,X_s) \mbox{d}W_s
- \end{equation}
- with
+ $$
+ with
\[
a(s,x) = \left( \begin{array}{@{\,}c@{\,}} 0 \\ g(s) \mu x_3 \\ \mu x_3 \end{array} \right),
\ \
@@ -494,7 +504,7 @@
solve.variable = c("x1", "x2", "x3"))
str(sabr.mod at parameter)
@
-The functions \code{f1}, \code{f2} and \code{f3} are defined in a way that, when the trajectory of the processes crosso zero from above, the trajectory is stopped ad zero. Notice that in this way the only visible parameter for \pkg{yuima} is \code{mu} as \code{rho} and \code{sig} are inside the bodies of the functions \code{f2} and \code{f3}. If we want to instruct \pkg{yuima} about these parameters, they should appear explicitly as arguments of the functions as explained by this \proglang{R} code
+The functions \code{f1}, \code{f2} and \code{f3} are defined in a way that, when the trajectory of the processes crosses zero from above, the trajectory is stopped at zero. Notice that in this way the only visible parameter for \pkg{yuima} is \code{mu} as \code{rho} and \code{sig} are inside the bodies of the functions \code{f2} and \code{f3}. If we want to instruct \pkg{yuima} about these parameters, they should appear explicitly as arguments of the functions as shown in the following \proglang{R} code
<<echo=TRUE, print=FALSE,results=hide>>=
f2 <- function(t, x1, x2, x3, rho, sig) {
ret <- 0
@@ -524,8 +534,8 @@
$$\de X_t = a(X_t) \de t + b(X_t) \de W_t^H
$$
where $W^H=\left( W^H_t, t\geq 0 \right)$ is a normalized fractional Brownian motion (fBM), {\it i.e.} the zero mean Gaussian processes with covariance function
-$$\dE( W_s^H W_t^H )= \frac{1}{2} \left[ |s|^{2H}+|t|^{2H}- |t-s|^{2H}\right]$$ with Hurst exponent $H\in(0,1)$. The fractional Brownian motion process is neither Markovian nor a semimartingale for $H\neq\frac{1}{2}$ but
-remains Gaussian. For $H > \frac12$, the solution $X_t$ above presents the long-range dependance property that makes it useful for different applications in biology, physics, ethernet traffic or finance.
+$$\dE( W_s^H W_t^H )= \frac{1}{2} \left( |s|^{2H}+|t|^{2H}- |t-s|^{2H}\right)$$ with Hurst exponent $H\in(0,1)$. The fractional Brownian motion process is neither Markovian nor a semimartingale for $H\neq\frac{1}{2}$ but
+remains Gaussian. For $H > \frac12$, the solution $X_t$ above presents the long-range dependence property that makes it useful for different applications in biology, physics, ethernet traffic or finance.
In order to specify a stochastic differential equation driven by fractional Gaussian noise it is necessary to specify the value of the Hurst parameter. For example, if we want to specify the following model
$$\de Y_t = 3 Y_t \de t + \de W_t^H$$
we proceed as follows
@@ -552,13 +562,29 @@
\subsection{L\'evy Processes}\label{sec:jump}
Jump processes can be specified in different ways in mathematics and hence in \pkg{yuima} package.
Let $Z_t$ be a Compound Poisson Process (i.e., jumps size follow some distribution, like the Gaussian law and jumps occur at Poisson times).
-Then it is possible to consider the following SDE which involves jumps
-$$\de X_t = a(X_t)\de t + b(X_t)\de W_t + \de Z_t$$
-In the next example we consider a compound Poisson process with intensity $\lambda=10$ with Gaussian jumps.
-This model can be specified in \code{setModel} using the argument \code{measure.type="CP"}
-A simple Ornstein-Uhlenbeck process with Gaussian jumps
+Then it is possible to consider the simple following SDE which involves jumps
+$$\de X_t = a(t,X_t,\theta)\de t + b(t,X_t,\theta)\de W_t + \de Z_t,$$
+or the more general representation
+\begin{equation}\label{eq:levy}
+\begin{aligned}
+\de X_t = & \,\,\, a(t,X_t,\theta)\de t + b(t,X_t,\theta)\de W_t + \int\limits_{|z|>1}\!\!\! c(X_{t-},z)\mu(\de t,\de z) \\
+&{}+\!\! \int\limits_{0<|z|\le 1}\!\!\! c(X_{t-},z)\{\mu(\de t,\de z) - \nu(\de z)\de t\}
+\end{aligned}, \quad X_0=x_0,
+\end{equation}
+with $\mu$ a random measure associated with the jumps of $X$, that is,
+$$
+\mu(dt,dz)=\sum_{s>0}\mathbf{1}_{\{\Delta Z_s \ne 0\}}\delta_{(s,\Delta Z_s)}(dt,dz),
+$$
+where $\delta$ denotes the Dirac measure and $Z$ the driving pure-jump L\'evy process of the form
+$$
+Z_{t}=\int_{0}^{t}\int_{|z|\le 1}z\{\mu(ds,dz) - \nu(dz)ds\}+\int_{0}^{t}\int_{|z|>1}z\mu(ds,dz).
+$$
+The L\'evy measure $\nu$ is any measure satisfying $\nu(\{0\})=0$ and $\int(1\wedge |z|^2)\nu(dz)< \infty$.
+In the next example we consider a compound Poisson process with intensity $\lambda=10$ with Gaussian jumps as driving L\'evy process.
+This model can be specified in \code{setModel} using the argument \code{measure.type="CP"} (fro Compound Poisson).
+A simple Ornstein-Uhlenbeck process with Gaussian jumps like this
$$\de X_t = -\theta X_t \de t + \sigma \de W_t + \de Z_t$$
-is specified as
+is then specified as follows
<<echo=TRUE, print=FALSE,fig=TRUE,width=9,height=4,results=hide>>=
mod5 <- setModel(drift=c("-theta*x"), diffusion="sigma",
jump.coeff="1", measure=list(intensity="10", df=list("dnorm(z, 0, 1)")),
@@ -569,7 +595,8 @@
@
\noindent
-Another possibility is to specify the jump component via the L\'evy measure. Without going into too much details, here is an example of specification of a simple Ornstein-Uhlenbeck process with IG (Inverse Gaussian) L\'evy measure
+Another possibility is to specify the jump component via the L\'evy measure $\nu(\cdot)$.
+ Without going into too much details, here is an example of specification of a simple Ornstein-Uhlenbeck process with IG (Inverse Gaussian) L\'evy measure $\nu(\cdot)$ and no Poisson component
$$\de X_t = -x \de t + \de Z_t$$
<<echo=TRUE, print=FALSE,fig=TRUE,width=9,height=4,results=hide>>=
mod6 <- setModel(drift="-x", xinit=1, jump.coeff="1",
@@ -579,25 +606,54 @@
plot(X)
@
-\subsection{Specification of Generic Models}
+\noindent
+The argument \code{type="code"} stands for user defined code.
+\subsection{Specification of Generic Models in \pkg{yuima}}
In general, the \pkg{yuima} package allows to specify a large family of models solutions to
-$$\de X_t =a(X_t)\de t + b(X_t)\de W_t + c(X_t)\de Z_t$$
+$$\de X_t =a(t,X_t,\theta)\de t + b(t,X_t,\theta)\de W_t + c(t,X_t,\theta)\de Z_t$$
using the following interface
<<echo=TRUE, print=FALSE,eval=FALSE>>=
setModel(drift, diffusion, hurst = 0.5, jump.coeff,
measure, measure.type, state.variable = "x",
jump.variable = "z", time.variable = "t", solve.variable, xinit)
@
+The coefficient $c()$ may also depend on the process $Z$, in this case the admissible form for the coefficient is $c(t,x,\theta, z) = c(t,x,\theta) \cdot z$. In the integral representation correspoding to formula \eqref{eq:levy} the coefficient $c(\cdot)$ is not allowed to depend on time though.
The \pkg{yuima} package implements many multivariate Random Numbers Generators (RNG) which are needed to simulate L\'evy paths including
\code{rIG} (Inverse Gaussian), \code{rNIG} (Normal Inverse Gaussian), \code{rbgamma} (Bilateral Gamma), \code{rngamma} (Normal Gamma) and \texttt{rstable} (Stable Laws).
Other user-defined RNG can be used freely.
\section{Simulation, Sampling and Subsampling}\label{sec:simul}
-The \code{simulate} function simulates \code{yuima} models according to Euler-Maruyama scheme in the presence of non-fractional diffusion noise and L\'evy jumps and uses the Cholesky or the Wood and Chan methods for the fractional Gaussian noise.
-The \code{simulate} function accepts several arguments including the description of the sampling structure, which is an object of type \code{yuima.sampling}. The \code{setSampling} allows for the specification of different sampling parameters including random sampling. Further, the \code{subsampling} allows us to subsample a trajectory of a simulated stochastic differential equation or a given time series in the \code{yuima.data} slot of a \code{yuima} object.
-Sampling and subsampling can be specified jointly as arguments to the \code{simulate} function. This is convenient if one wants to simulate data at very high frequency but then return only low frequency data for inference or other applications. In what follows we explain how to specify arguments of these \pkg{yuima} functions. Complete details can be found in the man pages of the \pkg{yuima} package.
+The \code{simulate} function simulates \code{yuima} models according to Euler-Maruyama scheme in the presence of non-fractional diffusion noise and L\'evy jumps and uses the Cholesky or the Wood and Chan methods for the fractional Gaussian noise. For diffusion models without jumps, \pkg{yuima} also implements the space discretized Euler-Maruyama method, which is more efficient, in the sense of mean squared error approximation, than the usual time-discretized Euler-Maruyama scheme.
+The (time discretized) Euler-Maruyama simulation scheme defines a grid of times
+$0 = \tau_0 < \tau_1 < \cdots < \tau_j < \tau_{j+1} < \cdots$,
+and constructs an approximative solution $\{\tilde{X}_t, t\geq 0\}$ of the original process $\{X_t, t\geq 0\}$ given by
+$$
+ \tilde{X}_{\tau_{j+1}} = \tilde{X}_{\tau_j} +
+ b(\tau_j,\tilde{X}_{\tau_j})(\tau_{j+1} - \tau_j) +
+ c(\tau_j, \tilde{X}_{\tau_j} )(W_{\tau_{j+1}}- W_{\tau_j}), \ \ j \geq 0.
+$$
+ As the interval $\tau_{j+1} - \tau_j$ is independent to $\{W_t\}_{t \geq \tau_j}$ for each $j$,
+ $W_{\tau_{j+1}} - W_{\tau_j}$ has the same distribution as $\sqrt{\tau_{j+1} - \tau_j}N_j$,
+ where $N_j$ is an i.i.d. sequence of standard normal variables.
+ The basic discretization scheme is equidistant, i.e.
+ $\tau_{j+1} - \tau_j = T/n = \Delta_n$, for all $j \geq 0$.
+But the \code{simulate} method provides also space-discretized Euler-Maruyama method to solve SDE.
+ This is at the moment designed for 1 factor SDE only, i.e., the case the driving
+ Brownian motion $W$ is 1-dimensional. In this case, the discretization scheme $\{\tau_j\}$
+ is defined as
+$$
+\tau_0 = 0, \ \ \tau_{j+1} = \inf\{ t > \tau_j; |W_t - W_{\tau_j}| = \epsilon \}
+$$ for each $j \geq 0$. Internally, \code{simulate} takes
+$\epsilon^2 = T/n = \Delta_n$
+which coincides with the mean of the interval $\tau_{j+1} - \tau_j$.
+This space-discretizing scheme is known to be 3 times efficient than the usual time-discretized scheme
+one in the sense of the mean squared error \citep{Fukasawa11}. To make use of the space-discretized Euler-Maruyama scheme, one should use the argument \code{space.discretized = TRUE} which by default, is set to \code{FALSE}.
+Other simulation scheme for 1-dimensional diffusion processes as explained in \citet{Iacus08} will be implemented in the near future.
-Assume that we want to simulate this model
+The \code{simulate} function accepts several arguments including the description of the sampling structure, which is an object of type \code{yuima.sampling}. The \code{setSampling} allows for the specification of different sampling parameters including random sampling. Further, the \code{subsampling} allows to subsample a trajectory of a simulated stochastic differential equation or a given time series in the \code{yuima.data} slot of a \code{yuima} object.
+Sampling and subsampling can be specified jointly as arguments to the \code{simulate} function. This is convenient if one wants to simulate data at very high frequency but then return only low frequency data for inference or other applications. In what follows we explain how to specify arguments of these \pkg{yuima} functions. Complete details can be found in the manual pages of the \pkg{yuima} package.
+
+Assume that we want to simulate the following model
$$
\begin{aligned}
\de X_{1,t} &= -\theta X_{1,t} \de t + \de W_{1,t} + X_{2,t} \de W_{3,t}\\
@@ -615,17 +671,17 @@
<<>>=
samp <- setSampling(Terminal=3, n=3000)
[TRUNCATED]
To get the complete diff run:
svnlook diff /svnroot/yuima -r 243
More information about the Yuima-commits
mailing list