R: Calculate the First Passage Time Density of a...
GQD.passage
R Documentation
Calculate the First Passage Time Density of a Time-Homogeneous GQD Process to a Fixed Barrier.
Description
GQD.passage uses the cumulant truncation procedure of Varughese (2013) in conjunction with a Volterra-type integral equation developed by Buonocore et al. (1987) in order to approximate the first passage time density of a time-homogeneous univariate GQD
Parameter vector giving the coefficients of the time-homogeneous GQD.
t
The time horizon up to and including which the density is to be evaluated.
delt
Stepsize for the solution of the first passage time density.
Value
density
The approximate first passage time density.
time
The time points at which the approximation is evaluated.
prob.coverage
The approximate cumulative probability coverage.
Warning
Warning [1]: Some instability may occur when delt is large or where the saddlepoint approximation fails. As allways it is important to check both the validity of the diffusion process under the given parameters and the quality of the sadlepoint approximation. For a given set of parameters the latter can be checked using GQD.density.
Warning [2]: The first passage time problem is considered from one side only i.e. Xs<B. For Xs>B one may equivalently consider Yt=-X_t, apply Ito's lemma and proceed as above.
Note
Note [1]: Since time-homogeneity is assumed for the present implementation, the interface of GQD.mcmc etc. is discarded and the model is inferred from the non-zero values of theta. For example theta = c(0.5*10,-0.5,0,0.5^2,0,0) defines an Ornstein-Uhlenbeck model.
A. Buonocore, A. Nobile, and L. Ricciardi. 1987 A new integral equation for the evaluation of first-passage-
time probability densities. Adv. Appl. Probab.19:784–800.
Daniels, H.E. 1954 Saddlepoint approximations in statistics. Ann. Math. Stat., 25:631–650.
Eddelbuettel, D. and Romain, F. 2011 Rcpp: Seamless R and C++ integration. Journal of Statistical Software, 40(8):1–18,. URL http://www.jstatsoft.org/v40/i08/.
Eddelbuettel, D. 2013 Seamless R and C++ Integration with Rcpp. New York: Springer. ISBN
978-1-4614-6867-7.
Eddelbuettel, D. and Sanderson, C. 2014 Rcpparmadillo: Accelerating r with high-performance C++
linear algebra. Computational Statistics and Data Analysis, 71:1054–1063. URL
http://dx.doi.org/10.1016/j.csda.2013.02.005.
Feagin, T. 2007 A tenth-order Runge-Kutta method with error estimate. In Proceedings of the IAENG
Conf. on Scientific Computing.
R. G. Jaimez, P. R. Roman and F. T. Ruiz. 1995 A note on the volterra integral equation for the
first-passage-time probability density. Journal of applied probability, 635–648.
Varughese, M.M. 2013 Parameter estimation for multivariate diffusion systems. Comput. Stat. Data An.,
57:417–428.
See Also
GQD.density for functions that generate the transitional density of GQDs. GQD.mcmc and GQD.remove for MCMC procedures to perform inference on GQDs.
Examples
#===============================================================================
# Calculate the first passage time from state X_0 = 7 to X_T =10 for
# various diffusions, with T the first passage time.
#===============================================================================
res1 <- GQD.passage(7,10,c(0.1*10,-0.1,0,0.2,0,0),25,1/100)
res2 <- GQD.passage(7,10,c(0.1*10,-0.1,0,0,0.2,0),25,1/100)
res3 <- GQD.passage(7,10,c(0,0.1*10,-0.1,0.5^2,0,0),25,1/100)
res4 <- GQD.passage(7,10,c(0,0.1*10,-0.1,0,0,0.1^2),25,1/100)
expr1 <- expression(dX[t]==(1-0.1*X[t])*dt+sqrt(0.2)*dW[t])
expr2 <- expression(dX[t]==(1-0.1*X[t])*dt+sqrt(0.2*X[t])*dW[t])
expr3 <- expression(dX[t]==(1*X[t]-0.1*X[t]^2)*dt+0.5*dW[t])
expr4 <- expression(dX[t]==(1*X[t]-0.1*X[t]^2)*dt+0.1*X[t]*dW[t])
#===============================================================================
# Plot the resulting first passage time densities.
#===============================================================================
par(mfrow=c(2,2))
plot(res1$density~res1$time,type='l',main=expr1,xlab='Time (t)',ylab='density',col='blue')
plot(res2$density~res2$time,type='l',main=expr2,xlab='Time (t)',ylab='density',col='blue')
plot(res3$density~res3$time,type='l',main=expr3,xlab='Time (t)',ylab='density',col='blue')
plot(res4$density~res4$time,type='l',main=expr4,xlab='Time (t)',ylab='density',col='blue')
Results
R version 3.3.1 (2016-06-21) -- "Bug in Your Hair"
Copyright (C) 2016 The R Foundation for Statistical Computing
Platform: x86_64-pc-linux-gnu (64-bit)
R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.
R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.
Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.
> library(DiffusionRgqd)
> png(filename="/home/ddbj/snapshot/RGM3/R_CC/result/DiffusionRgqd/GQD.passage.Rd_%03d_medium.png", width=480, height=480)
> ### Name: GQD.passage
> ### Title: Calculate the First Passage Time Density of a Time-Homogeneous
> ### GQD Process to a Fixed Barrier.
> ### Aliases: GQD.passage
> ### Keywords: syntax C++ first passage time
>
> ### ** Examples
>
> ## No test:
> #===============================================================================
> # Calculate the first passage time from state X_0 = 7 to X_T =10 for
> # various diffusions, with T the first passage time.
> #===============================================================================
>
> res1 <- GQD.passage(7,10,c(0.1*10,-0.1,0,0.2,0,0),25,1/100)
Compiling C++ code. Please wait. > res2 <- GQD.passage(7,10,c(0.1*10,-0.1,0,0,0.2,0),25,1/100)
Compiling C++ code. Please wait. > res3 <- GQD.passage(7,10,c(0,0.1*10,-0.1,0.5^2,0,0),25,1/100)
Compiling C++ code. Please wait. > res4 <- GQD.passage(7,10,c(0,0.1*10,-0.1,0,0,0.1^2),25,1/100)
Compiling C++ code. Please wait. >
> expr1 <- expression(dX[t]==(1-0.1*X[t])*dt+sqrt(0.2)*dW[t])
> expr2 <- expression(dX[t]==(1-0.1*X[t])*dt+sqrt(0.2*X[t])*dW[t])
> expr3 <- expression(dX[t]==(1*X[t]-0.1*X[t]^2)*dt+0.5*dW[t])
> expr4 <- expression(dX[t]==(1*X[t]-0.1*X[t]^2)*dt+0.1*X[t]*dW[t])
>
>
> #===============================================================================
> # Plot the resulting first passage time densities.
> #===============================================================================
>
> par(mfrow=c(2,2))
> plot(res1$density~res1$time,type='l',main=expr1,xlab='Time (t)',ylab='density',col='blue')
> plot(res2$density~res2$time,type='l',main=expr2,xlab='Time (t)',ylab='density',col='blue')
> plot(res3$density~res3$time,type='l',main=expr3,xlab='Time (t)',ylab='density',col='blue')
> plot(res4$density~res4$time,type='l',main=expr4,xlab='Time (t)',ylab='density',col='blue')
>
> ## End(No test)
>
>
>
>
>
> dev.off()
null device
1
>