Re: [stan-users] Question about Wiener Diffusion Model Distribution in rstan_2.7.0-2

[Bob Carpenter]

[Cross-posted to stan-dev]

Thanks for the detailed comments.

Could someone who understands this Wiener thing please respond
on stan-users?

I've never understood the function or its doc and we get a steady
stream of these messages saying it's wrong, which in the past
have been due to poor doc, not buggy code.

The next step would be to create reproducible errors.

- Bob

> On Feb 19, 2016, at 8:15 AM, Maciej Swat <macie...@gmail.com> wrote:
>
>
> Hi all
>
> I have a related comment.
> There are three potential typos in the STAN manual (2.9.0) in the
> PDF for 'Wiener Diffusion Model Distributions’, but I’m not sure if these
> bugs 2&3 proliferated into the code.
>
> 1. y is defined as 'y \in (0, \tau)’ but it should be 'y \in (\tau, \infty)’
> otherwise the term \sqrt{y - \tau} is undefined.
>
>
> Two others based on reference paper (Blurton et at. 2012)
>
> 2. the exponent of the denominator (y - \tau) should be 3/2 instead of 3
> See attached PDF plots for
> alphaB=3.25;tauB=0.25;betaB=0.4;deltaB=0.65; N = 150;
>
> 3. the numerator in standard normal distribution function \phi reads
> ‘(2k + \beta) \alpha’ instead of '2k + \alpha \beta’.
> The difference between these two results in very minor changes.
>
> Best Regards, Maciej
>
>
> --
> You received this message because you are subscribed to the Google Groups "Stan users mailing list" group.
> To unsubscribe from this group and stop receiving emails from it, send an email to stan-users+...@googlegroups.com.
> To post to this group, send email to stan-...@googlegroups.com.
> For more options, visit https://groups.google.com/d/optout.
> <Rplot03.pdf><Rplot05.pdf>
>
>
> --
> You received this message because you are subscribed to the Google Groups "Stan users mailing list" group.
> To unsubscribe from this group and stop receiving emails from it, send an email to stan-users+...@googlegroups.com.
> To post to this group, send email to stan-...@googlegroups.com.
> For more options, visit https://groups.google.com/d/optout.
> <Screen Shot 2016-02-19 at 10.44.08.png>
>
> --
> You received this message because you are subscribed to the Google Groups "Stan users mailing list" group.
> To unsubscribe from this group and stop receiving emails from it, send an email to stan-users+...@googlegroups.com.
> To post to this group, send email to stan-...@googlegroups.com.
> For more options, visit https://groups.google.com/d/optout.
> <Screen Shot 2016-02-19 at 10.44.19.png>
>
>> On 28 Jan 2016, at 20:30, Royce Anders <ander...@gmail.com> wrote:
>>
>> ah right.. I think the trick is in what I missed before, the 1-b, with the -v, and the Wiener dist pdf here is only calculating the probability that the upper bound is hit, as acting like the lower bound due to inversing these two parameters together. Great :)
>>
>> On Thu, Jan 28, 2016 at 6:30 PM, Mike Lawrence <mike....@gmail.com> wrote:
>> My recollection is that the distribution is defined so that negative values represent the times of opposite-boundary crossings. So you code your data with negative signs on the rts for error trials, and the distribution knows to treat these as opposite boundary crossings.
>>
>> Mike
>>
>>
>> On Thu, Jan 28, 2016 at 1:05 PM, Royce Anders <ander...@gmail.com> wrote:
>> in this implementation of the two-boundary Wiener process in STAN/JAGS, is it missing a piece of information: particularly when the process terminates due to hitting the opposite boundary than the sign of the drift rate?
>>
>> e.g., : if (resp[n] == 1) {
>> y[n] ~ wiener(a, t, b, v);
>> }
>> else {
>> y[n] ~ wiener(a, t, 1-b, -v);
>> }
>>
>> Is it not a completely true two-boundary Wiener process then?
>>
>> I understand that the probability is quite low that the opposite boundary (e.g., the lower boundary will be hit when the drift v is positive), however it becomes increasingly probable as the threshold values (a) are small, or the drifts (v) are slow. In the case of (a) it may be quite important to capture speed-accuracy manipulations in experiments. It's my understanding that this is an important/defining aspect of the diffusion decision model (DDM), hence am I missing something here or is this mirrored application (e.g., as in two conjugate inverse gaussians) capturing this important aspect of the data when we fit by this implementation?
>>
>>
>> On Friday, July 31, 2015 at 1:57:10 PM UTC+2, Andrew Ellis wrote:
>> this seems to work:
>>
>> lower boundary responses should no longer be coded as negative response times; instead the response times and responses should be passed to stan as data:
>>
>> library(rstan)
>> library(RWiener)
>>
>> set.seed(0)
>> dat <- rwiener(n = 5000, alpha = 1, tau = 0.2, beta = 0.75, delta = 0.5)
>> standata <- list(y = dat$q,
>> resp = ifelse(dat$resp == "upper", 1, 0),
>> N = length(dat$resp))
>>
>>
>> and in the stan model, y ~ wiener(a, t, b, v) if the "upper" response was selected, and y ~ wiener(a, t, 1-b, -v) if the "lower" response was selected.
>>
>> data {
>> int<lower=0> N;
>> real y[N];
>> int<lower=0, upper=1> resp[N];
>> }
>> parameters {
>> real<lower=0> a;
>> real<lower=0> t;
>> real<lower=0, upper=1> b;
>> real v;
>> }
>> model {
>> v ~ normal(0, 2);
>> a ~ gamma(3, 3);
>> t ~ gamma(4, 10);
>> b ~ beta(5, 5);
>>
>>
>> for (n in 1:N) {
>>
>>
>> if (resp[n] == 1) {
>> y[n] ~ wiener(a, t, b, v);
>> }
>> else {
>> y[n] ~ wiener(a, t, 1-b, -v);
>> }
>>
>>
>> }
>> }
>>
>> Full code example:
>>
>> library(rstan)
>> library(RWiener)
>>
>>
>> set.seed(0)
>> dat <- rwiener(n = 5000, alpha = 1, tau = 0.2, beta = 0.75, delta = 0.5)
>>
>> ## stan model ----
>> wiener_model <- "
>> data {
>> int<lower=0> N;
>> real y[N];
>> int<lower=0, upper=1> resp[N];
>> }
>> parameters {
>> real<lower=0> a;
>> real<lower=0> t;
>> real<lower=0, upper=1> b;
>> real v;
>> }
>> model {
>> v ~ normal(0, 2);
>> a ~ gamma(3, 3);
>> t ~ gamma(4, 10);
>> b ~ beta(5, 5);
>>
>>
>> for (n in 1:N) {
>>
>>
>> if (resp[n] == 1) {
>> y[n] ~ wiener(a, t, b, v);
>> }
>> else {
>> y[n] ~ wiener(a, t, 1-b, -v);
>> }
>>
>>
>> }
>> }"
>>
>>
>> ## stan data ----
>> standata <- list(y = dat$q,
>> resp = ifelse(dat$resp == "upper", 1, 0),
>> N = length(dat$resp))
>>
>> fit <- stan(model_code = wiener_model, data = standata)
>>
>>
>>
>>
>> On Thursday, July 30, 2015 at 12:50:12 PM UTC+2, Andrew Ellis wrote:
>> Has anyone successfully used this?
>>
>> I tried to implement the example given here: http://www.cidlab.com/supp.php?o=wienerstan
>>
>> The 'wiener.hpp' file states that "y should contain reaction times in seconds, with upper-boundary responses strictly positive and lower-boundary response times coded as strictly negative numbers."
>>
>> My example code:
>>
>> library(rstan)
>> library(RWiener)
>>
>>
>> # generate data using RWiener package
>> set.seed(0)
>> dat <- rwiener(n = 1000, alpha = 1, tau = 0.2, beta = 0.5, delta = 0.5)
>>
>> standata <- list(y = ifelse(dat$resp == "upper", dat$q, -1 * dat$q),
>> N = length(dat$resp))
>>
>>
>> wiener_model <- "
>> data {
>> int<lower=0> N;
>> real y[N];
>> }
>> parameters {
>> real<lower=0> a;
>> real<lower=0> t;
>> real<lower=0,upper=1> b;
>> real v;
>> }
>> model {
>> v ~ normal(0, 2);
>>
>>
>> for (n in 1:N) {
>> y[n] ~ wiener(1, 0.2, 0.5, v);
>> }
>> }"
>>
>>
>> fit <- stan(model_code = wiener_model, data = standata)
>>
>>
>>
>> I get the error message:
>> Rejecting initial value:
>> Error evaluating the log probability at the initial value.
>>
>>
>> Informational Message: The current Metropolis proposal is about to be rejected because of the following issue:
>> Exception thrown at line 16:
>> stan::math::wiener_log(%1%): Random variable is -0.395369, but must be > 0!
>> If this warning occurs sporadically, such as for highly constrained variable types like covariance matrices, then the sampler is fine,
>> but if this warning occurs often then your model may be either severely ill-conditioned or misspecified.
>>
>>
>>
>> According to the manual (p. 396), Stan returns the first passage time of the accumulation process over the upper boundary only. How do I deal with both decision boundaries, if I cannot pass in lower-boundary response times as negative numbers?
>>
>> cheers,
>> Andrew
>>
>>
>>
>>> sessionInfo()
>> R version 3.2.1 (2015-06-18)
>> Platform: x86_64-apple-darwin13.4.0 (64-bit)
>> Running under: OS X 10.10.4 (Yosemite)
>>
>>
>> locale:
>> [1] en_GB.UTF-8/en_GB.UTF-8/en_GB.UTF-8/C/en_GB.UTF-8/en_GB.UTF-8
>>
>>
>> attached base packages:
>> [1] stats graphics grDevices utils datasets methods base
>>
>>
>> other attached packages:
>> [1] RWiener_1.2-0 rstan_2.7.0-2 Rcpp_0.12.0 magrittr_1.5 devtools_1.8.0
>>
>>
>> loaded via a namespace (and not attached):
>> [1] codetools_0.2-14 digest_0.6.8 R6_2.1.0 stats4_3.2.1 git2r_0.10.1 httr_1.0.0
>> [7] stringi_0.5-5 curl_0.9.1 xml2_0.1.1 tools_3.2.1 stringr_1.0.0.9000 inline_0.3.14
>> [13] rversions_1.0.2 memoise_0.2.1
>>
>>
>> --
>> You received this message because you are subscribed to the Google Groups "Stan users mailing list" group.
>> To unsubscribe from this group and stop receiving emails from it, send an email to stan-users+...@googlegroups.com.
>> To post to this group, send email to stan-...@googlegroups.com.
>> For more options, visit https://groups.google.com/d/optout.
>>
>>
>> --
>> You received this message because you are subscribed to a topic in the Google Groups "Stan users mailing list" group.
>> To unsubscribe from this topic, visit https://groups.google.com/d/topic/stan-users/-6wJfA-t2cQ/unsubscribe.
>> To unsubscribe from this group and all its topics, send an email to stan-users+...@googlegroups.com.
>> To post to this group, send email to stan-...@googlegroups.com.
>> For more options, visit https://groups.google.com/d/optout.
>>
>>
>> --
>> You received this message because you are subscribed to the Google Groups "Stan users mailing list" group.
>> To unsubscribe from this group and stop receiving emails from it, send an email to stan-users+...@googlegroups.com.
>> To post to this group, send email to stan-...@googlegroups.com.
>> For more options, visit https://groups.google.com/d/optout.
>
> --
> You received this message because you are subscribed to the Google Groups "Stan users mailing list" group.
> To unsubscribe from this group and stop receiving emails from it, send an email to stan-users+...@googlegroups.com.
> To post to this group, send email to stan-...@googlegroups.com.
> For more options, visit https://groups.google.com/d/optout.