Monte-Carlo Simulation and Kernel Density Estimation of First passage time

A.C. Guidoum1 and K. Boukhetala2

2020-03-27

The fptsdekd() functions

A new algorithm based on the Monte Carlo technique to generate the random variable FPT of a time homogeneous diffusion process (1, 2 and 3D) through a time-dependent boundary, order to estimate her probability density function.

Let \(X_t\) be a diffusion process which is the unique solution of the following stochastic differential equation:

\[\begin{equation}\label{eds01} dX_t = \mu(t,X_t) dt + \sigma(t,X_t) dW_t,\quad X_{t_{0}}=x_{0} \end{equation}\]

if \(S(t)\) is a time-dependent boundary, we are interested in generating the first passage time (FPT) of the diffusion process through this boundary that is we will study the following random variable:

\[ \tau_{S(t)}= \left\{ \begin{array}{ll} inf \left\{t: X_{t} \geq S(t)|X_{t_{0}}=x_{0} \right\} & \hbox{if} \quad x_{0} \leq S(t_{0}) \\ inf \left\{t: X_{t} \leq S(t)|X_{t_{0}}=x_{0} \right\} & \hbox{if} \quad x_{0} \geq S(t_{0}) \end{array} \right. \]

The main arguments to ‘random’ fptsdekd() (where k=1,2,3) consist:

The following statistical measures (S3 method) for class fptsdekd() can be approximated for F.P.T \(\tau_{S(t)}\):

The main arguments to ‘density’ dfptsdekd() (where k=1,2,3) consist:

Examples

FPT for 1-Dim SDE

Consider the following SDE and linear boundary:

\[\begin{align*} dX_{t}= & (1-0.5 X_{t}) dt + dW_{t},~x_{0} =1.7.\\ S(t)= & 2(1-sinh(0.5t)) \end{align*}\]

Generating the first passage time (FPT) of this model through this boundary: \[ \tau_{S(t)}= \inf \left\{t: X_{t} \geq S(t) |X_{t_{0}}=x_{0} \right\} ~~ \text{if} \quad x_{0} \leq S(t_{0}) \]

Set the model \(X_t\):

Generate the first-passage-time \(\tau_{S(t)}\), with fptsde1d() function ( based on density() function in [base] package):

Itô Sde 1D:
    | dX(t) = (1 - 0.5 * X(t)) * dt + 1 * dW(t)
    | t in [0,1].
Boundary:
    | S(t) = 2 * (1 - sinh(0.5 * t))
F.P.T:
    | T(S(t),X(t)) = inf{t >=  0 : X(t) >=  2 * (1 - sinh(0.5 * t)) }
    | Crossing realized 965 among 1000.
 [1] 0.107542 0.203707 0.089800 0.067849 0.240718 0.187920 0.416042
 [8] 0.054746 0.025669 0.061076

The following statistical measures (S3 method) for class fptsde1d() can be approximated for the first-passage-time \(\tau_{S(t)}\):

[1] 0.20638
[1] 0.041248
[1] 0.13043
[1] 0.061989
       0%       25%       50%       75%      100% 
0.0097554 0.0603962 0.1304262 0.2780265 0.9962872 
[1] 5.0625
[1] 1.5834
[1] 0.98461
[1] 0.0097554
[1] 0.99629
[1] 0.0086314
[1] 0.031954

The result summaries of the first-passage-time \(\tau_{S(t)}\):


Monte-Carlo Statistics of F.P.T:
|T(S(t),X(t)) = inf{t >=  0 : X(t) >=  2 * (1 - sinh(0.5 * t)) }
                        
Mean             0.20638
Variance         0.04129
Median           0.13043
Mode             0.06199
First quartile   0.06040
Third quartile   0.27803
Minimum          0.00976
Maximum          0.99629
Skewness         1.58338
Kurtosis         5.06251
Coef-variation   0.98461
3th-order moment 0.01329
4th-order moment 0.00863
5th-order moment 0.00499
6th-order moment 0.00320

Display the exact first-passage-time \(\tau_{S(t)}\), see Figure 1:

The kernel density approximation of ‘fpt1d’, using dfptsde1d() function (hist=TRUE based on truehist() function in MASS package), see e.g. Figure 2.

Since fptdApprox and DiffusionRgqd packages can very effectively handle first passage time problems for diffusions with analytically tractable transitional densities we use it to compare some of the results from the Sim.DiffProc package.

fptsde1d() vs Approx.fpt.density()

Consider for example a diffusion process with SDE:

\[\begin{align*} dX_{t}= & 0.48 X_{t} dt + 0.07 X_{t} dW_{t},~x_{0} =1.\\ S(t)= & 7 + 3.2 t + 1.4 t \sin(1.75 t) \end{align*}\] The resulting object is then used by the Approx.fpt.density() function in package fptdApprox to approximate the first passage time density:

Using fptsde1d() and dfptsde1d() functions in the Sim.DiffProc package:

Itô Sde 1D:
    | dX(t) = 0.48 * X(t) * dt + 0.07 * X(t) * dW(t)
    | t in [0,10].
Boundary:
    | S(t) = 7 + 3.2 * t + 1.4 * t * sin(1.75 * t)
F.P.T:
    | T(S(t),X(t)) = inf{t >=  0 : X(t) >=  7 + 3.2 * t + 1.4 * t * sin(1.75 * t) }
    | Crossing realized 1000 among 1000.
 [1] 8.4117 5.8765 6.2605 6.0622 6.0488 8.3695 6.0661 6.4495 6.1632
[10] 8.4051

Monte-Carlo Statistics of F.P.T:
|T(S(t),X(t)) = inf{t >=  0 : X(t) >=  7 + 3.2 * t + 1.4 * t * sin(1.75 * t) }
                         
Mean              6.51683
Variance          0.88330
Median            6.11985
Mode              6.06244
First quartile    5.96986
Third quartile    6.37719
Minimum           5.49777
Maximum           9.10482
Skewness          1.49664
Kurtosis          3.51957
Coef-variation    0.14422
3th-order moment  1.24246
4th-order moment  2.74605
5th-order moment  5.30548
6th-order moment 10.83898

By plotting the approximations:

 `fptsde1d()` vs `Approx.fpt.density()`

fptsde1d() vs Approx.fpt.density()

fptsde1d() vs GQD.TIpassage()

Consider for example a diffusion process with SDE:

\[\begin{align*} dX_{t}= & \theta_{1}X_{t}(10+0.2\sin(2\pi t)+0.3\sqrt(t)(1+\cos(3\pi t))-X_{t}) ) dt + \sqrt(0.1) X_{t} dW_{t},~x_{0} =8.\\ S(t)= & 12 \end{align*}\] The resulting object is then used by the GQD.TIpassage() function in package DiffusionRgqd to approximate the first passage time density:

Using fptsde1d() and dfptsde1d() functions in the Sim.DiffProc package:

Itô Sde 1D:
    | dX(t) = theta1 * X(t) * (10 + 0.2 * sin(2 * pi * t) + 0.3 * sqrt(t) *     (1 + cos(3 * pi * t)) - X(t)) * dt + sqrt(0.1) * X(t) * dW(t)
    | t in [1,4].
Boundary:
    | S(t) = 12
F.P.T:
    | T(S(t),X(t)) = inf{t >=  1 : X(t) >=  12 }
    | Crossing realized 924 among 1000.
 [1] 1.2159 1.4396 2.9346 3.3582 2.9095 2.0872 2.1483 2.6700 2.1994
[10] 1.5246

Monte-Carlo Statistics of F.P.T:
|T(S(t),X(t)) = inf{t >=  1 : X(t) >=  12 }
                        
Mean             2.14945
Variance         0.48577
Median           2.06170
Mode             1.46280
First quartile   1.51230
Third quartile   2.65939
Minimum          1.11332
Maximum          3.97248
Skewness         0.56845
Kurtosis         2.27088
Coef-variation   0.32426
3th-order moment 0.19246
4th-order moment 0.53586
5th-order moment 0.48835
6th-order moment 0.88108

By plotting the approximations (hist=TRUE based on truehist() function in MASS package):

`fptsde1d()` vs `GQD.TIpassage()`

fptsde1d() vs GQD.TIpassage()

FPT for 2-Dim SDE’s

The following \(2\)-dimensional SDE’s with a vector of drift and a diagonal matrix of diffusion coefficients:

\[\begin{equation}\label{eq:09} \begin{cases} dX_t = f_{x}(t,X_{t},Y_{t}) dt + g_{x}(t,X_{t},Y_{t}) dW_{1,t}\\ dY_t = f_{y}(t,X_{t},Y_{t}) dt + g_{y}(t,X_{t},Y_{t}) dW_{2,t} \end{cases} \end{equation}\]

\(W_{1,t}\) and \(W_{2,t}\) is a two independent standard Wiener process. First passage time (2D) \((\tau_{(S(t),X_{t})},\tau_{(S(t),Y_{t})})\) is defined as:

\[ \left\{ \begin{array}{ll} \tau_{S(t),X_{t}}=\inf \left\{t: X_{t} \geq S(t)|X_{t_{0}}=x_{0} \right\} & \hbox{if} \quad x_{0} \leq S(t_{0}) \\ \tau_{S(t),Y_{t}}=\inf \left\{t: Y_{t} \geq S(t)|Y_{t_{0}}=y_{0} \right\} & \hbox{if} \quad y_{0} \leq S(t_{0}) \end{array} \right. \] and \[ \left\{ \begin{array}{ll} \tau_{S(t),X_{t}}= \inf \left\{t: X_{t} \leq S(t)|X_{t_{0}}=x_{0} \right\} & \hbox{if} \quad x_{0} \geq S(t_{0}) \\ \tau_{S(t),Y_{t}}= \inf \left\{t: Y_{t} \leq S(t)|Y_{t_{0}}=y_{0} \right\} & \hbox{if} \quad y_{0} \geq S(t_{0}) \end{array} \right. \]

Assume that we want to describe the following Stratonovich SDE’s (2D):

\[\begin{equation}\label{eq016} \begin{cases} dX_t = 5 (-1-Y_{t}) X_{t} dt + 0.5 Y_{t} \circ dW_{1,t}\\ dY_t = 5 (-1-X_{t}) Y_{t} dt + 0.5 X_{t} \circ dW_{2,t} \end{cases} \end{equation}\]

and \[ S(t)=\sin(2\pi t) \]

Set the system \((X_t , Y_t)\):

Generate the couple \((\tau_{(S(t),X_{t})},\tau_{(S(t),Y_{t})})\), with fptsde2d() function::

Stratonovich Sde 2D:
    | dX(t) = 5 * (-1 - Y(t)) * X(t) * dt + 0.5 * Y(t) o dW1(t)
    | dY(t) = 5 * (-1 - X(t)) * Y(t) * dt + 0.5 * X(t) o dW2(t)
    | t in [0,1].
Boundary:
    | S(t) = sin(2 * pi * t)
F.P.T:
    | T(S(t),X(t)) = inf{t >=  0 : X(t) <=  sin(2 * pi * t) }
    |   And 
    | T(S(t),Y(t)) = inf{t >=  0 : Y(t) >=  sin(2 * pi * t) }
    | Crossing realized 1000 among 1000.
         x       y
1  0.11048 0.51071
2  0.12736 0.50116
3  0.12843 0.50947
4  0.13097 0.50479
5  0.12776 0.50534
6  0.12599 0.51098
7  0.15425 0.50238
8  0.13745 0.49875
9  0.15475 0.50511
10 0.12286 0.49946

The following statistical measures (S3 method) for class fptsde2d() can be approximated for the couple \((\tau_{(S(t),X_{t})},\tau_{(S(t),Y_{t})})\):

[1] 0.13326 0.50353
[1] 0.000181883 0.000028843
[1] 0.13276 0.50337
[1] 0.13302 0.50240
$x
      0%      25%      50%      75%     100% 
0.091453 0.124124 0.132761 0.141578 0.181096 

$y
     0%     25%     50%     75%    100% 
0.48677 0.50007 0.50337 0.50711 0.52155 
[1] 3.2311 3.1736
[1]  0.324641 -0.024634
[1] 0.101253 0.010671
[1] 0.091453 0.486773
[1] 0.18110 0.52155
[1] 0.0000001071029 0.0000000026455
[1] 0.00033529 0.06432632

The result summaries of the couple \((\tau_{(S(t),X_{t})},\tau_{(S(t),Y_{t})})\):


Monte-Carlo Statistics for the F.P.T of (X(t),Y(t))
    | T(S(t),X(t)) = inf{t >=  0 : X(t) <=  sin(2 * pi * t) }
    |    And
    | T(S(t),Y(t)) = inf{t >=  0 : Y(t) >=  sin(2 * pi * t) }
                  T(S,X)   T(S,Y)
Mean             0.13326  0.50353
Variance         0.00018  0.00003
Median           0.13276  0.50337
Mode             0.13302  0.50240
First quartile   0.12412  0.50007
Third quartile   0.14158  0.50711
Minimum          0.09145  0.48677
Maximum          0.18110  0.52155
Skewness         0.32464 -0.02463
Kurtosis         3.23110  3.17359
Coef-variation   0.10125  0.01067
3th-order moment 0.00000  0.00000
4th-order moment 0.00000  0.00000
5th-order moment 0.00000  0.00000
6th-order moment 0.00000  0.00000

Display the exact first-passage-time \(\tau_{S(t)}\), see Figure 5:

The marginal density of \((\tau_{(S(t),X_{t})}\) and \(\tau_{(S(t),Y_{t})})\) are reported using dfptsde2d() function, see e.g. Figure 6.

A contour and image plot of density obtained from a realization of system \((\tau_{(S(t),X_{t})},\tau_{(S(t),Y_{t})})\).

A \(3\)D plot of the Joint density with:

Return to fptsde2d()

FPT for 3-Dim SDE’s

The following \(3\)-dimensional SDE’s with a vector of drift and a diagonal matrix of diffusion coefficients:

\[\begin{equation}\label{eq17} \begin{cases} dX_t = f_{x}(t,X_{t},Y_{t},Z_{t}) dt + g_{x}(t,X_{t},Y_{t},Z_{t}) dW_{1,t}\\ dY_t = f_{y}(t,X_{t},Y_{t},Z_{t}) dt + g_{y}(t,X_{t},Y_{t},Z_{t}) dW_{2,t}\\ dZ_t = f_{z}(t,X_{t},Y_{t},Z_{t}) dt + g_{z}(t,X_{t},Y_{t},Z_{t}) dW_{3,t} \end{cases} \end{equation}\] \(W_{1,t}\), \(W_{2,t}\) and \(W_{3,t}\) is a 3 independent standard Wiener process. First passage time (3D) \((\tau_{(S(t),X_{t})},\tau_{(S(t),Y_{t})},\tau_{(S(t),Z_{t})})\) is defined as:

\[ \left\{ \begin{array}{ll} \tau_{S(t),X_{t}}=\inf \left\{t: X_{t} \geq S(t)|X_{t_{0}}=x_{0} \right\} & \hbox{if} \quad x_{0} \leq S(t_{0}) \\ \tau_{S(t),Y_{t}}=\inf \left\{t: Y_{t} \geq S(t)|Y_{t_{0}}=y_{0} \right\} & \hbox{if} \quad y_{0} \leq S(t_{0}) \\ \tau_{S(t),Z_{t}}=\inf \left\{t: Z_{t} \geq S(t)|Z_{t_{0}}=z_{0} \right\} & \hbox{if} \quad z_{0} \leq S(t_{0}) \end{array} \right. \] and \[ \left\{ \begin{array}{ll} \tau_{S(t),X_{t}}= \inf \left\{t: X_{t} \leq S(t)|X_{t_{0}}=x_{0} \right\} & \hbox{if} \quad x_{0} \geq S(t_{0}) \\ \tau_{S(t),Y_{t}}= \inf \left\{t: Y_{t} \leq S(t)|Y_{t_{0}}=y_{0} \right\} & \hbox{if} \quad y_{0} \geq S(t_{0}) \\ \tau_{S(t),Z_{t}}= \inf \left\{t: Z_{t} \leq S(t)|Z_{t_{0}}=z_{0} \right\} & \hbox{if} \quad z_{0} \geq S(t_{0}) \\ \end{array} \right. \]

Assume that we want to describe the following SDE’s (3D): \[\begin{equation}\label{eq0166} \begin{cases} dX_t = 4 (-1-X_{t}) Y_{t} dt + 0.2 dW_{1,t}\\ dY_t = 4 (1-Y_{t}) X_{t} dt + 0.2 dW_{2,t}\\ dZ_t = 4 (1-Z_{t}) Y_{t} dt + 0.2 dW_{3,t} \end{cases} \end{equation}\] and \[ S(t)=-1.5+3t \]

Set the system \((X_t , Y_t , Z_t)\):

Generate the triplet \((\tau_{(S(t),X_{t})},\tau_{(S(t),Y_{t})},\tau_{(S(t),Z_{t})})\), with fptsde3d() function::

Itô Sde 3D:
    | dX(t) = 4 * (-1 - X(t)) * Y(t) * dt + 0.2 * dW1(t)
    | dY(t) = 4 * (1 - Y(t)) * X(t) * dt + 0.2 * dW2(t)
    | dZ(t) = 4 * (1 - Z(t)) * Y(t) * dt + 0.2 * dW3(t)
    | t in [0,1].
Boundary:
    | S(t) = -1.5 + 3 * t
F.P.T:
    | T(S(t),X(t)) = inf{t >=  0 : X(t) <=  -1.5 + 3 * t }
    |   And 
    | T(S(t),Y(t)) = inf{t >=  0 : Y(t) >=  -1.5 + 3 * t }
    |   And 
    | T(S(t),Z(t)) = inf{t >=  0 : Z(t) <=  -1.5 + 3 * t }
    | Crossing realized 1000 among 1000.
         x        y       z
1  0.52670 0.024297 0.78568
2  0.50634 0.023111 0.80254
3  0.54032 0.023006 0.78718
4  0.52996 0.024950 0.76237
5  0.52370 0.023213 0.80861
6  0.55294 0.023462 0.80313
7  0.54789 0.023400 0.75685
8  0.53502 0.024272 0.80064
9  0.54227 0.023340 0.75972
10 0.52042 0.026061 0.74419

The following statistical measures (S3 method) for class fptsde3d() can be approximated for the triplet \((\tau_{(S(t),X_{t})},\tau_{(S(t),Y_{t})},\tau_{(S(t),Z_{t})})\):

[1] 0.53251 0.02320 0.78277
[1] 0.0001832797 0.0000016427 0.0008785594
[1] 0.53202 0.02316 0.78254
[1] 0.529438 0.023076 0.781245
$x
     0%     25%     50%     75%    100% 
0.49597 0.52343 0.53202 0.54113 0.57705 

$y
      0%      25%      50%      75%     100% 
0.019666 0.022292 0.023160 0.024105 0.027957 

$z
     0%     25%     50%     75%    100% 
0.68368 0.76365 0.78254 0.80254 0.87476 
[1] 2.9722 2.9424 3.2112
[1]  0.186492  0.096446 -0.191374
[1] 0.025436 0.055273 0.037885
[1] 0.495971 0.019666 0.683675
[1] 0.577047 0.027957 0.874757
[1] 0.0000001000405227 0.0000000000079564 0.0000024836149789
[1] 0.08072072765 0.00000029503 0.37865036624

The result summaries of the triplet \((\tau_{(S(t),X_{t})},\tau_{(S(t),Y_{t})},\tau_{(S(t),Z_{t})})\):


Monte-Carlo Statistics for the F.P.T of (X(t),Y(t),Z(t))
    | T(S(t),X(t)) = inf{t >=  0 : X(t) <=  -1.5 + 3 * t }
    |    And
    | T(S(t),Y(t)) = inf{t >=  0 : Y(t) >=  -1.5 + 3 * t }
    |    And
    | T(S(t),Z(t)) = inf{t >=  0 : Z(t) <=  -1.5 + 3 * t }
                  T(S,X)  T(S,Y)   T(S,Z)
Mean             0.53251 0.02320  0.78277
Variance         0.00018 0.00000  0.00088
Median           0.53202 0.02316  0.78254
Mode             0.52944 0.02308  0.78125
First quartile   0.52343 0.02229  0.76365
Third quartile   0.54113 0.02410  0.80254
Minimum          0.49597 0.01967  0.68368
Maximum          0.57705 0.02796  0.87476
Skewness         0.18649 0.09645 -0.19137
Kurtosis         2.97220 2.94245  3.21124
Coef-variation   0.02544 0.05527  0.03789
3th-order moment 0.00000 0.00000  0.00000
4th-order moment 0.00000 0.00000  0.00000
5th-order moment 0.00000 0.00000  0.00000
6th-order moment 0.00000 0.00000  0.00000

Display the exact first-passage-time \(\tau_{S(t)}\), see Figure 9:

The marginal density of \(\tau_{(S(t),X_{t})}\) ,\(\tau_{(S(t),Y_{t})}\) and \(\tau_{(S(t),Z_{t})})\) are reported using dfptsde3d() function, see e.g. Figure 10.


Marginal density for the F.P.T of X(t)
    | T(S,X) = inf{t >= 0 : X(t) <= -1.5 + 3 * t}

Data: out[, "x"] (1000 obs.);   Bandwidth 'bw' = 0.0029862

       x                f(x)        
 Min.   :0.48701   Min.   : 0.0015  
 1st Qu.:0.51176   1st Qu.: 0.9375  
 Median :0.53651   Median : 5.7959  
 Mean   :0.53651   Mean   :10.0917  
 3rd Qu.:0.56126   3rd Qu.:19.2485  
 Max.   :0.58601   Max.   :29.0495  

Marginal density for the F.P.T of Y(t)
    | T(S,Y) = inf{t >= 0 : Y(t) >= -1.5 + 3 * t}

Data: out[, "y"] (1000 obs.);   Bandwidth 'bw' = 0.0002899

       y                 f(y)        
 Min.   :0.018796   Min.   :  0.015  
 1st Qu.:0.021304   1st Qu.:  4.939  
 Median :0.023812   Median : 45.848  
 Mean   :0.023812   Mean   : 99.598  
 3rd Qu.:0.026319   3rd Qu.:201.593  
 Max.   :0.028827   Max.   :311.998  

Marginal density for the F.P.T of Z(t)
    | T(S,Z) = inf{t >= 0 : Z(t) <= -1.5 + 3 * t}

Data: out[, "z"] (1000 obs.);   Bandwidth 'bw' = 0.0065603

       z                 f(z)        
 Min.   :0.018796   Min.   :  0.015  
 1st Qu.:0.021304   1st Qu.:  4.939  
 Median :0.023812   Median : 45.848  
 Mean   :0.023812   Mean   : 99.598  
 3rd Qu.:0.026319   3rd Qu.:201.593  
 Max.   :0.028827   Max.   :311.998  

For an approximate joint density for \((\tau_{(S(t),X_{t})},\tau_{(S(t),Y_{t})},\tau_{(S(t),Z_{t})})\) (for more details, see package sm or ks.)

Return to fptsde3d()

Further reading

  1. snssdekd() & dsdekd() & rsdekd()- Monte-Carlo Simulation and Analysis of Stochastic Differential Equations.
  2. bridgesdekd() & dsdekd() & rsdekd() - Constructs and Analysis of Bridges Stochastic Differential Equations.
  3. fptsdekd() & dfptsdekd() - Monte-Carlo Simulation and Kernel Density Estimation of First passage time.
  4. MCM.sde() & MEM.sde() - Parallel Monte-Carlo and Moment Equations for SDEs.
  5. TEX.sde() - Converting Sim.DiffProc Objects to LaTeX.
  6. fitsde() - Parametric Estimation of 1-D Stochastic Differential Equation.

References

  1. Boukhetala K (1996). Modelling and Simulation of a Dispersion Pollutant with Attractive Centre, volume 3, pp. 245-252. Computer Methods and Water Resources, Computational Mechanics Publications, Boston, USA.

  2. Boukhetala K (1998). Estimation of the first passage time distribution for a simulated diffusion process. Maghreb Mathematical Review, 7, pp. 1-25.

  3. Boukhetala K (1998). Kernel density of the exit time in a simulated diffusion. The Annals of The Engineer Maghrebian, 12, pp. 587-589.

  4. Guidoum AC, Boukhetala K (2020). Sim.DiffProc: Simulation of Diffusion Processes. R package version 4.5, URL https://cran.r-project.org/package=Sim.DiffProc.

  5. Pienaar EAD, Varughese MM (2016). DiffusionRgqd: An R Package for Performing Inference and Analysis on Time-Inhomogeneous Quadratic Diffusion Processes. R package version 0.1.3, URL https://CRAN.R-project.org/package=DiffusionRgqd.

  6. Roman, R.P., Serrano, J. J., Torres, F. (2008). First-passage-time location function: Application to determine first-passage-time densities in diffusion processes. Computational Statistics and Data Analysis. 52, 4132-4146.

  7. Roman, R.P., Serrano, J. J., Torres, F. (2012). An R package for an efficient approximation of first-passage-time densities for diffusion processes based on the FPTL function. Applied Mathematics and Computation, 218, 8408-8428.


  1. Department of Probabilities & Statistics, Faculty of Mathematics, University of Science and Technology Houari Boumediene, BP 32 El-Alia, U.S.T.H.B, Algeria, E-mail ()

  2. Faculty of Mathematics, University of Science and Technology Houari Boumediene, BP 32 El-Alia, U.S.T.H.B, Algeria, E-mail ()