Description Usage Arguments Details Value Warning Note Author(s) References See Also Examples
GQD.TIpassage()
solves first passage time problems for GQDs with time dependent coefficients:
dX_t = (G_0(t)+G_1(t)X_t+G_2(t)X_t^2)dt+√{Q_0(t)+Q_1(t)X_t+Q_2(t)X_t^2}dW_t
to a fixed barrier. The function combines the cumulant truncation procedure of Varughese (2013) with a numerical solution to a non-singular Volterra integral equation for the first passage time density, developed by Buonocore et al. (1987), in order to generate an approximate solution.
1 | GQD.TIpassage(Xs,B, s, t, delt, theta=c(0), IEQ.type='Buonocore', wrt=FALSE)
|
Xs |
Initial value of the diffusion process at time |
B |
Fixed barrier (or first constant in static barier transform - see detail [1]). |
s |
Starting time for the diffusion process (see detail [2]). |
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. |
theta |
Parameter vector for parameters contained in the coefficients (if required). |
IEQ.type |
Currently only |
wrt |
If |
Detail [1]: First passage throug a time dependant barrier may be analised by applying the transform:
Y_t = X_t -B_t,
if B_t may can be decomposed as
B_t = k+f(t).
By applying Ito's lemma the revised drift and diffusion functionals, and first passage parameters may be inferred.
Detail [2]: The starting time is of particular importance when the drift and/or diffusion terms are time-inhomogeneous. Take care to select the correct starting time - especially if the drift or diffusion components whch are time dependant have poles or singular points in the time domain.
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 [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 saddlepoint 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 [1]: The coefficients od the GQD may be parameterized using the reserved variable theta
. For example:
G0 <- function(t){theta[1]*(theta[2]+sin(2*pi*(t-theta[3])))}
may be used so long as values are asigned in the function call, say
GQD.TIpassage(Xs=3,B=11,tmin=1,tmax=10,delt=1/100,theta=c(1,10,0.5))
.
Note [2]: Due to syntactical differences between R and C++ special functions have to be used when terms that depend on t
. When the function cannot be separated in to terms that contain a single t
, the prod(a,b)
function must be used. For example (see examples below):
G0 <- function(t){0.1*(10+0.2*sin(2*pi*t)+0.3*prod(sqrt(t),1+cos(3*pi*t)))}
.
Here sqrt(t)*cos(3*pi*t) constitutes the product of two terms that cannot be written i.t.o. a single t
. To circumvent this isue, one may use the prod(a,b)
function.
Note [3]: Similarly, the ^ - operator is not overloaded in C++. Instead the pow(x,p)
function may be used to calculate x^p. For example sin(2*pi*t)^3 in:
G0 <- function(t){0.1*(10+0.2*pow(sin(2*pi*t),3))}
.
Note [4]: delt
is used as the stepsize of the Runge-Kutta method used to numerically solve a system of ODEs used to approximate the cumulants of the underlying diffusion procees. The 10th-order scheme of Feagin (2007) is used as the default method.
Etienne A.D. Pienaar: etiennead@gmail.com
Updates available on GitHub at https://github.com/eta21.
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.
GQD.density
for functions that generate the transitional density of GQDs and GQD.mcmc
for MCMC procedures for GQD models.
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 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 | #=========================================================================
# Generate the first passage time density of a CIR process with time
# dependant drift to a fixed barrier.
#-------------------------------------------------------------------------
# Remove any existing coefficients.
GQD.remove()
# Define the coefficients of the process.
G0 <- function(t){10+0.5*sin(2*pi*t)}
G1 <- function(t){-1}
Q1 <- function(t){0.25}
#Define the parameters of the first passage time problem.
delt <- 1/100 # The stepsize for the solution
X0 <- 8 # The initial value of the process
BB <- 11 # Fixed barrier
T0 <- 2 # Starting time of the diffusion
TT <- 10 # Time horizon of the computation
# Run the calculation
res1 <- GQD.TIpassage(X0,BB,T0,TT,delt)
# Remove any existing coefficients.
GQD.remove()
# Redefine the coefficients.
G0 <- function(t){ 0.1*(10+0.2*sin(2*pi*t)+0.3*prod(sqrt(t),1+cos(3*pi*t)))}
G1 <- function(t){-0.1*(1+0.2*sin(2*pi*t))}
Q1 <- function(t){0.25}
# Redefine the parameters of the f.p.t. problem.
delt <- 1/100
X0 <- 8
BB <- 11
T0 <- 1
TT <- 10
# Run the calculation
res2 <- GQD.TIpassage(X0,BB,T0,TT,delt)
#===============================================================================
# Plot the two solutions.
#===============================================================================
expr1 <- expression(dX[t]==(10+0.5*sin(2*pi*t)-X[t])*dt+0.25*sqrt(X[t])*dW[t])
expr2 <- expression(dX[t]==(0.1*(10+0.2*sin(2*pi*t)+0.3*sqrt(t)*(1+cos(3*pi*t))
-0.1*(1+0.2*sin(2*pi*t))*X[t])*dt+0.25*sqrt(X[t])*dW[t]))
par(mfrow=c(1,1))
plot(res1$density~res1$time,type='l',col='blue',
ylab='Density',xlab='Time',main=expr1,cex.main=0.95)
plot(res2$density~res2$time,type='l',col='red',
ylab='Density',xlab='Time',main =expr2,cex.main=0.95)
#===============================================================================
# Let's see how sensitive the first passage density is w.r.t a speed parameter
# of a non-linear diffusion.
#===============================================================================
GQD.remove()
# Redefine the coefficients with a parameter theta:
G1 <- function(t){theta[1]*(10+0.2*sin(2*pi*t)+0.3*prod(sqrt(t),1+cos(3*pi*t)))}
G2 <- function(t){-theta[1]}
Q2 <- function(t){0.1}
# Now just give a value for the parameter in the standard fashion:
res3=GQD.TIpassage(8,12,1,4,1/100,theta=c(0.5))
plot(res3$density~res3$time,type='l',col=2,ylim=c(0,1.0),
main='First Passage Time Density',ylab='Density',xlab='Time',cex.main=0.95)
# Change the parameter and see the effect on the f.p.t. density.
th.seq=seq(0.1,0.5,1/20)
for(i in 2:length(th.seq))
{
res3=GQD.TIpassage(8,12,1,4,1/100,,theta=c(th.seq[i]))
lines(res3$density~res3$time,type='l',col=rainbow(10)[i])
}
lines(res3$density~res3$time,type='l',col=rainbow(10)[i],lwd=2)
legend('topright',legend=th.seq,col=rainbow(10),lty='solid',cex=0.75,
title=expression(theta[1]))
#===============================================================================
|
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.