Documentation for Wiener first passage time distribution

Modeling
1

(This question may be more to do with some minor issues with the Stan documentation, and so I presumably should open an issue on the documentation’s GitHub repository. But before I do so, I think it is probably better to ask a question here.)

On the documentation page for the Wiener first passage time distribution, the density function of first passage times is stated as

\text{Wiener}(y|\alpha, \tau, \beta, \delta) = \frac{\alpha^3}{(y-\tau)^{3/2}} \exp \! \left(- \delta \alpha \beta - \frac{\delta^2(y-\tau)}{2}\right) \sum_{k = - \infty}^{\infty} (2k + \beta) \phi \! \left(\frac{2k \alpha + \beta}{\sqrt{y - \tau}}\right),

where \phi is the standard normal density.

Given that Joachim Vandekerckhove implemented the C++ code for Stan’s wiener functions (see math/wiener_lpdf.hpp at 3a080d74172c19f7eb9e7339214bcc4e2847c3bb · stan-dev/math · GitHub), I presume https://journal.r-project.org/archive/2014-1/vandekerckhove-wabersich.pdf is a definitive reference. The above function appears to be close to, but not exactly, what you obtain when we substitute \frac{t-\tau}{\alpha^2} into Equation 2 (in Appendix), and then substitute the result back as the last term of Equation 1, and substitute y for t. When I do that, I get the following:

\frac{1}{\sqrt{2\pi}} \frac{\alpha}{(y-\tau)^{3/2}}\exp \! \left(- \alpha \beta \delta - \frac{1}{2} \delta^2(y-\tau)\right) \sum_{k=-\infty}^{\infty} (2k + \beta) \exp{\left(-\frac{1}{2} \frac{\alpha^2 (\beta + 2k)^2}{y-\tau}\right)},

and writing that using the standard normal density, I get

\frac{\alpha}{(y-\tau)^{3/2}}\exp \! \left(- \alpha \beta \delta - \frac{1}{2} \delta^2(y-\tau)\right) \sum_{k=-\infty}^{\infty} (2k + \beta) \phi{\left(\frac{\alpha (\beta + 2k)}{\sqrt{y-\tau}}\right)}.

That’s close to the formula above, but not identical. In particular, we have \alpha rather than \alpha^2 as the numerator of the first term, and we have \alpha(\beta + 2k) rather than 2k\alpha +\beta as the numerator inside \phi.

Have I made an embarrassing error with my calculations or is the error in the density formula in the Stan documentation?

On a related note, Eq 1 in https://journal.r-project.org/archive/2014-1/vandekerckhove-wabersich.pdf is the bivariate density function over first passage times over the lower barrier in a two barrier Wiener diffusion process. In the C++ code, it is definitely implementing the first passage time for the upper barrier. It is implementing the method described by Navarro and Fuss (2009), which is gives the first passage for the lower barrier, but they replace \beta with 1-\beta and \delta for -\delta to give the first passage time of the upper barrier. So we know, by looking at the C++, what the Stan code is implementing, but is the formula in the documentation not giving (notwithstanding any errors) the density for the lower barrier? (It is also clearly not the case that the density in the documentation is Eq 1 of Wabersich and Vandekerckhove with \beta and \delta exchanged for 1-\beta and -\delta, respectively.)

As a final point, would it not be best to make it explicit, as Wabersich and Vandekerckhove do, that the density function is the joint density of first passage time and the identity of the barrier. Wabersich and Vandekerckhove write the lhs of their formula as such: d(t, x = 0 \vert \alpha, \tau, \beta, \delta), making it clear it is a bivariate density. As it is in Stan documentation, it appears to be presented as a univariate density, and so one whose integral from y=\tau to y=\infty is 1.0, which it is not.

2 Likes
2

Hi @mark-andrews, apologies for this being left unanswered. I’m going through Wiener first passage issues and noticed it. Did you ever create an issue or resolve the discrepancy?

1 Like