[Yuima-commits] r250 - in pkg: yuima/R yuimadocs/inst/doc/JSS

noreply at r-forge.r-project.org noreply at r-forge.r-project.org
Thu Oct 10 17:00:25 CEST 2013


Author: iacus
Date: 2013-10-10 17:00:23 +0200 (Thu, 10 Oct 2013)
New Revision: 250

Added:
   pkg/yuimadocs/inst/doc/JSS/final-edited.Rnw
Modified:
   pkg/yuima/R/mmfrac.R
   pkg/yuima/R/qgv.R
Log:
print methods

Modified: pkg/yuima/R/mmfrac.R
===================================================================
--- pkg/yuima/R/mmfrac.R	2013-10-10 14:54:34 UTC (rev 249)
+++ pkg/yuima/R/mmfrac.R	2013-10-10 15:00:23 UTC (rev 250)
@@ -74,7 +74,10 @@
 }
 
 print.mmfrac<-function(x,...){
-x
+tmp <- rbind(x$coefficients, diag(x$vcov))
+rownames(tmp) <- c("Estimate", "Std. Error")
+cat("\nFractional OU estimation\n")
+print(tmp)
 }
 
 

Modified: pkg/yuima/R/qgv.R
===================================================================
--- pkg/yuima/R/qgv.R	2013-10-10 14:54:34 UTC (rev 249)
+++ pkg/yuima/R/qgv.R	2013-10-10 15:00:23 UTC (rev 250)
@@ -222,6 +222,11 @@
 
 }
 
+
 print.qgv<-function(x,...){
-x
+    tmp <- rbind(x$coefficients, diag(x$vcov))
+    rownames(tmp) <- c("Estimate", "Std. Error")
+    cat("\nFractional OU estimation\n")
+    print(tmp)
 }
+

Added: pkg/yuimadocs/inst/doc/JSS/final-edited.Rnw
===================================================================
--- pkg/yuimadocs/inst/doc/JSS/final-edited.Rnw	                        (rev 0)
+++ pkg/yuimadocs/inst/doc/JSS/final-edited.Rnw	2013-10-10 15:00:23 UTC (rev 250)
@@ -0,0 +1,1752 @@
+\documentclass[article]{jss}
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%% declarations for jss.cls %%
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+%\usepackage{Sweave}    
+%\usepackage{Rd}    
+
+
+%%% INFOS FOR YUIMA CORE TEAM %%%
+
+% USE THIS IF YOU DON'T WANT TO EXECUTE R CODE
+\SweaveOpts{prefix.string=yuima, echo=TRUE, eval=FALSE}
+
+% USE THIS INSTEAD IF YOU WANT TO EXECUTE R CODE
+\SweaveOpts{prefix.string=yuima, echo=TRUE, eval=TRUE}
+
+%% before editing this file get the new version, type this command on Terminal
+%% in the same directory where this Rnw-file lives
+%% svn up
+
+%% HOW TO BUILD THIS FILE
+%% from inside the directory where this Rnw-file lives
+%% ./script.sh article-new
+
+%% to commit to svn the changes type
+%% from inside the directory where this Rnw-file lives
+%% svn commit -m "my changes"
+
+\usepackage{amsmath,amssymb}
+
+\usepackage[utf8x]{inputenc} 
+
+\usepackage{comment}
+\excludecomment{en-text}
+
+\usepackage{enumerate}%\begin{enumerate}[{(}i{)}]
+\usepackage{bm}
+\usepackage{bbm} %\mathbbm{ABCabc\nabla}
+\usepackage[bbgreekl]{mathbbol} %\mathbb{abcABC\nabla\bbalpha\bbbeta\Delta  }
+\usepackage{natbib}
+\usepackage{dsfont}
+\usepackage{amssymb}
+\usepackage{amsthm}
+\usepackage{latexsym}
+\usepackage{amsmath}
+
+\def\be{\begin{equation}}
+\def\ee{\end{equation}}
+\def\bea{\begin{eqnarray}}
+\def\eea{\end{eqnarray}}
+\def\beas{\begin{eqnarray*}}
+\def\eeas{\end{eqnarray*}}
+\def\sskip{\hspace{5mm}}
+
+\newcommand{\colorr}{\color[rgb]{0.8,0,0}}
+\newcommand{\colorg}{\color[rgb]{0,0.5,0}}
+\newcommand{\colorb}{\color[rgb]{0,0,0.8}}
+\newcommand{\colord}{\color[rgb]{0.8,0.3,0}}
+\newcommand{\colorn}{\color[rgb]{1,1,1}}
+\newcommand{\colorred}{\color[rgb]{0.8,0,0}}
+\newcommand{\colory}{\color{yellow}}
+\newcommand{\coloro}{\color[rgb]{1,0.4,0}}%{1,0.851,0}}
+\newcommand{\coloroy}{\color[rgb]{1,0.95,0}}
+\newcommand{\colorsb}{\color[rgb]{0,0.95,1}}
+\newcommand{\colorro}{\color[rgb]{0.851,0.255,0.467}} %rose
+
+
+%% almost as usual
+\author{Alexandre Brouste\\University of Le Mans \And 
+        Masaaki Fukasawa\\Osaka 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
+        Hiroki Masuda\\Kyushu University \And
+        Ryosuke Nomura\\University of Tokyo\AND
+        Teppei Ogihara\\Osaka University \And
+        Yasutaka Shimuzu\\Osaka University \And
+        Masayuki Uchida\\Osaka University \And
+        Nakahiro Yoshida\\University of Tokyo 
+        }
+\title{The YUIMA Project: a Computational Framework for Simulation and Inference of Stochastic Differential Equations}
+
+%% for pretty printing and a nice hypersummary also set:
+\Plainauthor{Alexandre Brouste, Stefano M. Iacus, Hiroki Masuda, Masayuki Uchida, Nakahiro Yoshida} %% comma-separated
+\Plaintitle{YUIMA: Simulation and Inference for SDE} %% without formatting
+\Shorttitle{The YUIMA Project} %% a short title (if necessary)
+
+%% an abstract and keywords
+\Abstract{
+The Yuima Project is an open source and collaborative effort  aimed at developing the \proglang{R} package named \pkg{yuima} for simulation and inference of stochastic differential equations. 
+
+In the \pkg{yuima} package stochastic differential equations can be of very abstract type, multidimensional, driven by Wiener process or fractional Brownian motion with general Hurst parameter, with or without jumps specified as L\'evy noise. 
+
+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 provides some examples of applications.
+}
+\Keywords{inference for stochastic processes, simulation, stochastic differential equations}
+\Plainkeywords{inference for stochastic processes, simulation, stochastic differential equations} %% without formatting
+%% at least one keyword must be supplied
+
+%% publication information
+%% NOTE: Typically, this can be left commented and will be filled out by the technical editor
+%% \Volume{13}
+%% \Issue{9}
+%% \Month{September}
+%% \Year{2004}
+%% \Submitdate{2004-09-29}
+%% \Acceptdate{2004-09-29}
+
+%% The address of (at least) one author should be given
+%% in the following format:
+\Address{
+  Stefano M. Iacus\\
+  Department of Economics, Management and Quantitative Methods\\
+  University of Milan\\
+  Via Conservatorio 7, 20122 Milan, Italy\\
+  E-mail: \email{stefano.iacus at unimi.it}\\
+  URL: \url{http://www.economia.unimi.it/~iacus/}
+\\
+\\
+  Nakahiro Yoshida\\
+  Graduate School of Mathematical Science\\
+  University of Tokyo\\
+  3-8-1 Komaba, Meguro-ku,  Tokyo 153-8914, Japan\\
+  E-mail: \email{nakahiro at ms.u-tokyo.ac.jp}\\
+  URL: \url{http://www.ms.u-tokyo.ac.jp/~nakahiro/hp-naka-e}
+}
+%% It is also possible to add a telephone and fax number
+%% before the e-mail in the following format:
+%% Telephone: +43/1/31336-5053
+%% Fax: +43/1/31336-734
+
+
+%% preamble for Rnw files
+<<print=FALSE, echo=FALSE>>=
+options(prompt="R> ")
+options(width=60)
+@
+
+\def\ep{{\varepsilon}}
+\def\ve{{\varepsilon}}
+\def\de{{\rm d}}
+\def\dE{{\mathbb E}}
+\def\diag{\mathop{\rm diag}\nolimits}
+
+%% for those who use Sweave please include the following line (with % symbols):
+%% need no \usepackage{Sweave.sty}
+%% end of declarations %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+
+\begin{document}
+
+
+
+
+%% include your article here, just as usual
+%% Note that you should use the \pkg{}, \proglang{} and \code{} commands.
+
+%\section[About Java]{About \proglang{Java}}
+%% Note: If there is markup in \(sub)section, then it has to be escape as above.
+
+
+
+
+
+
+
+
+\section{Introduction}
+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
+to carry out simulation and statistical analysis. 
+Both categories of procedures may depend on each other.
+Statistical inference often requires a simulation technique 
+as a subroutine, and a certain simulation method 
+needs to fix a tuning parameter by applying 
+a statistical methodology. 
+It is especially the case of stochastic processes 
+because most of expected values involved 
+do not admit an explicit expression. 
+The \pkg{yuima} package facilitates comprehensive, systematic approaches 
+to the solution. 
+
+
+Stochastic differential equations are 
+commonly used 
+to model random evolution along continuous or 
+practically continuous time, such as 
+the random movements of stock prices. 
+Theory of statistical inference for 
+stochastic differential equations already has a fairly long history, 
+more than three decades, but it is still developing quickly new 
+methodologies and expanding the area. 
+The formulas produced by the theory are usually very sophisticated, 
+which makes it difficult for standard users not necessarily 
+familiar with this field to enjoy utilities. 
+For example, by taking advantage of the analytic approach,  the asymptotic expansion method for computing 
+option prices (i.e., expectation of an irregular functional of 
+a stochastic process) provides precise approximation values 
+instantaneously. The expansion formula, which has a long expression involving more than 900 terms of multiple integrals, is already coded in the \pkg{yuima} package for generic diffusion processes.
+
+
+The \pkg{yuima} package delivers up-to-date methods as a package 
+onto the desk of the user working 
+with simulation and/or statistics for stochastic differential equations. 
+In the \pkg{yuima} package, stochastic differential equations can be of very abstract type, 
+multidimensional, driven by Wiener process or fractional Brownian motion 
+with general Hurst parameter, with or without jumps specified as L\'evy noise. 
+
+
+
+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{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.
+
+
+Although we assume that the reader of this paper has a basic knowledge of the \proglang{R} language, most of the examples are easy to understand 
+if he/she knows stochastic differential equations intuitively or symbolically. 
+
+
+\section{The \pkg{yuima} Package}\label{sec2}
+The package \pkg{yuima} depends on some other packages, like  \pkg{zoo} \citep{zoo}, which can be installed separately.
+The package \pkg{zoo} is used internally to store time series data. This dependence may change in the future adopting a more flexible class for internal storage of time series.
+
+\subsection{How to Obtain the Package}
+The \pkg{yuima} package is hosted on R-Forge and the web page
+of the Project is  \url{http://r-forge.r-project.org/projects/yuima}.
+The  R-Forge page contains the latest development version, and stable
+version of the package will also be available through CRAN.
+Development versions of the package are not supposed to be stable or functional, thus
+the occasional user should consider to install the stable version first.
+The package can be installed from R-Forge using
+\code{install.packages("yuima",repos="http://R-Forge.R-project.org")}.
+If, the R-Forge system does not provide binary builds of the \pkg{yuima} package, the
+user can also try 
+\code{install.packages("yuima",repos="http://R-Forge.R-project.org", type="source")}.
+
+
+\subsection{The Main Object and Classes}
+Before discussing the methods for simulation and inference for stochastic processes solutions to stochastic differential equations, here we discuss the main classes in the package.
+As mentioned there are different classes of objects defined in the \pkg{yuima} package and
+the main class is called the \code{yuima-class}. This class  is composed of several slots.
+Figure \ref{fig:classes} represents the different classes and their slots.
+\begin{figure}[th]
+\centering{\includegraphics[width=\textwidth]{yuima-class}}
+\caption{The main classes in the \pkg{yuima} package.}
+\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 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 (by ``\code{yuima} object'' we mean a, possibly incomplete, object of class \code{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 (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 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
+$$ \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
+$$
+\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{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 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. Further, \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. 
+All functions in the package are assumed to get as much information
+as possible from the model instead of replicating the same code everywhere.
+If there are missing pieces of information, we may change or extend the description
+of the model.
+\\
+\\
+An object of class \code{yuima.model} contains several slots listed below. 
+To see inside its structure, one can use the \proglang{R} command \code{str} on a \code{yuima} object.
+\begin{itemize}
+\item \code{drift} is an \proglang{R} vector of expressions which contains the drift specification;
+\item \code{diffusion} is itself a list of 1 slot which describes the diffusion
+  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 (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;
+\item \code{common} contains the names of the parameters in common between the drift and diffusion coefficients;
+\item \code{diffusion} contains the parameters belonging to the diffusion coefficient;
+\item \code{drift} contains the parameters belonging to the drift coefficient;
+\item \code{jump} contains the parameters belonging to the  coefficient of the L\'evy noise;
+\item \code{measure} contains the parameters describing the L\'evy measure (explained details in Section \ref{sec:jump});
+\end{itemize}
+\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 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 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} 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} 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} 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. (independent and identically distributed) sampling while estimators of parameters in the drift coefficient are slower or, in some cases, not even consistent. The \pkg{yuima} package 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 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
+$$\de X_t = -3 X_t \de t + \frac{1}{1+X_t^2}\de W_t.$$
+In the above $a(x) = -3 x$ and $b(x) = \frac{1}{1+x^2}$ according to the notation of previous section and $W_t$ is a standard Wiener process.
+This can be described in \pkg{yuima} by specifying the drift and diffusion coefficients as plain
+\proglang{R} expressions passed as strings
+<<print=FALSE,echo=FALSE,results=hide>>=
+library(yuima)
+@ 
+<<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 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 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)
+plot(X)
+@
+
+\noindent
+The \code{simulate} function fills in addition the two slots \code{data} and \code{sampling} of the \code{yuima} object. 
+<<>>=
+str(X,vec.len=2)
+@
+More details on how to change the default sampling scheme for the \code{simulate} method and how to perform subsampling will be given in Section \ref{sec:simul}.
+
+
+\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 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
+<<echo=TRUE, print=FALSE,results=hide>>=
+mod1b <- setModel(drift = "-3*s*y", diffusion = "1/(1+y^2)", 
+state.var="y", time.var="s")
+@
+In this case the solution variable is the same as the state variable. Indeed, the \code{yuima.model} object appears as follows
+<<>>=
+str(mod1b)
+@
+Once again, the user may use the \code{simulate} method to perform simulation.
+
+\subsection{Specification of Parametric Models}\label{sec:par}
+Assume now that we want to specify a parametric model like this
+$$\de X_t = -\theta X_t \de t + \frac{1}{1+X_t^\gamma}\de W_t$$
+where $a(x,\theta) = -\theta x$ and $b(x,\gamma) = 1/(1+x^\gamma)$.
+Then, \pkg{yuima} attempts to distinguish the parameters' names from the ones of the state and time variables
+<<echo=TRUE, print=FALSE,results=hide>>=
+mod2 <- setModel(drift = "-theta*x", diffusion = "1/(1+x^gamma)")
+@
+so, in this case, \code{theta} and \code{gamma}, which are different form \code{x} and \code{t}, are assumed to be parameters. Notice that in the above notation $\theta$ and $\gamma$ are generic names for the components of a parameters' vector $\theta$ in the notation of Section \ref{sec:model}.
+<<>>=
+str(mod2)
+@
+In order to simulate the parametric model it is necessary to specify the values of the parameters $\theta$ and $\gamma$ as shown in the next code chunk
+<<echo=TRUE, print=FALSE,fig=TRUE,height=4,results=hide>>=
+set.seed(123)
+X <- simulate(mod2,true.param=list(theta=1,gamma=3))
+plot(X)
+@
+
+\subsection{Multidimensional Processes}\label{sec:multi}
+Next is an example of a system of two stochastic differential equations for the couple $(X_{1,t}, X_{2,t})$ driven by three independent Brownian motions $(W_{1,t}, W_{2,t}, W_{3,t})$
+$$
+\begin{aligned}
+\de X_{1,t} &= -3 X_{1,t} \de t + \de W_{1,t} + X_{2,t} \de W_{3,t}\\
+\de X_{2,t} &= -(X_{1,t} + 2 X_{2,t}) \de t + X_{1,t} \de W_{1,t} + 3 \de W_{2,t}
+\end{aligned}
+$$
+but this has to be organized into matrix form with a vector of drift expressions and a diffusion matrix
+$$
+\left(\begin{array}{c}\de X_{1,t} \\\de X_{2,t}\end{array}\right)=
+\left(\begin{array}{c} -3 X_{1,t} \\  -X_{1,t} - 2X_{2,t}\end{array}\right)\de t +
+\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 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 
+b <- matrix(c("1","x1","0","3","x2","0"),2,3)  #  diffusion matrix
+mod3 <- setModel(drift = a, diffusion = b, solve.variable = sol)
+@
+Looking at the structure of the \code{noise.number} slot in \code{mod3}, one can see that this is now set to 3 which is taken as the number of columns of the diffusion matrix.
+Again, this model can be easily simulated
+<<echo=TRUE, print=FALSE,fig=TRUE,width=9,height=4,results=hide>>=
+set.seed(123)
+X <- simulate(mod3)
+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{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}
+$$
+$$ (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
+ $$
+ 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
+ $$
+ with
+ \[
+  a(s,x) = \left( \begin{array}{@{\,}c@{\,}} 0 \\ g(s) \mu x_3 \\ \mu x_3 \end{array} \right),
+ \ \ 
+ b(s,x) = \left( \begin{array}{@{\,}cc@{\,}} x_2 | x_1|^{2/3} &  0 \\
+  g(s)\sigma \rho x_3 &  g(s)\sigma \sqrt{1-\rho^2} x_3 \\ 
+ \sigma \rho  x_3 & \sigma \sqrt{1-\rho^2} x_3 
+ \end{array} \right)
+ \]
+ for $ x = (x_1, x_2, x_3)$. 
+<<echo=TRUE, print=FALSE,results=hide>>=
+mu <- 0.1
+sig <- 0.2
+rho <- -0.7
+
+g <- function(t) {0.4 + (0.1 + 0.2*t)* exp(-2*t)}
+
+
+f1 <- function(t, x1, x2, x3) {
+    ret <- 0
+    if(x1 > 0 && x2 > 0) ret <- x2*exp(log(x1)*2/3)
+    return(ret)
+}
+
+f2 <- function(t, x1, x2, x3) {
+     ret <- 0
+     if(x3 > 0) ret <- rho*sig*x3
+     return(ret)
+}
+
+f3 <- function(t, x1, x2, x3) {
+     ret <- 0
+     if(x3 > 0) ret <- sqrt(1-rho^2)*sig*x3
+     return(ret)
+}
+
+diff.coef.matrix <- matrix(c("f1(t,x1,x2,x3)",
+ "f2(t,x1,x2,x3) * g(t)", "f2(t,x1,x2,x3)", "0", 
+     "f3(t,x1,x2,x3)*g(t)", "f3(t,x1,x2,x3)"),  3, 2)
+
+sabr.mod <- setModel(drift = c("0", "mu*g(t)*x3", "mu*x3"), 
+diffusion = diff.coef.matrix, state.variable = c("x1", "x2", "x3"),
+    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 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
+     if(x3 > 0) ret <- rho*sig*x3
+     return(ret)
+ }
+
+ f3 <- function(t, x1, x2, x3, rho, sig) {
+     ret <- 0
+     if(x3 > 0) ret <- sqrt(1-rho^2)*sig*x3
+     return(ret)
+ }
+
+ diff.coef.matrix <- matrix(c("f1(t,x1,x2,x3)", 
+ "f2(t,x1,x2,x3,rho, sig) * g(t)", "f2(t,x1,x2,x3,rho,sig)", 
+ "0", "f3(t,x1,x2,x3,rho,sig)*g(t)", "f3(t,x1,x2,x3,rho,sig)"),  3, 2)
+
+ sabr.mod <- setModel(drift = c("0", "mu*g(t)*x3", "mu*x3"), 
+ diffusion = diff.coef.matrix, state.variable = c("x1", "x2", "x3"), 
+ solve.variable = c("x1", "x2", "x3"))
+str(sabr.mod at parameter)
+@
+
+
+\subsection{Fractional Gaussian Noise}\label{sec:fgn}
+The \pkg{yuima} allows for the description of stochastic differential equations driven by fractional Brownian motion of the following type
+$$\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 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
+<<echo=TRUE, print=FALSE,fig=TRUE,width=9,height=6,results=hide>>=
+mod4A <- setModel(drift="3*y", diffusion=1, hurst=0.3, solve.var="y")
+mod4B <- setModel(drift="3*y", diffusion=1, hurst=0.7, solve.var="y")
+set.seed(123)
+X1 <- simulate(mod4A,sampling=setSampling(n=1000))
+X2 <- simulate(mod4B,sampling=setSampling(n=1000))
+par(mfrow=c(2,1))
+par(mar=c(2,3,1,1))
+plot(X1,main="H=0.3")
+plot(X2,main="H=0.7")
+@
+
+\noindent
+In this case, the appropriate slot is now filled
+@ 
+<<>>=
[TRUNCATED]

To get the complete diff run:
    svnlook diff /svnroot/yuima -r 250


More information about the Yuima-commits mailing list