itsdfish / SequentialSamplingModels.jl

A unified interface for simulating and evaluating sequential sampling models in Julia.
https://itsdfish.github.io/SequentialSamplingModels.jl/dev/
MIT License
27 stars 4 forks source link

Negative RT DDM implementation and misc #21

Closed itsdfish closed 1 year ago

itsdfish commented 1 year ago

Hi @kiante-fernandez. I am moving a discussion of the PR here. Unfortunately, I couldn't figure out how to view the PR locally before merging it (it's possible that it would have worked, but I was in the wrong folder).

I would like to bring up some issues. First, the constructors were not working. I ended up commenting out the following to get it to work:

# function DDM(ν::T, α::T, τ::T, z::T; check_args::Bool=true) where {T <: Real}
#     check_args && Distributions.@check_args DDM (α, α > zero(α)) (τ, τ > zero(τ)) (z, z ≥ 0 && z ≤ 1)
#     return DDM{T}(ν, α, τ, z)
# end

# function DDM(ν::T, α::T, τ::T; check_args::Bool=true) where {T <: Real}
#     return DDM(ν, α, τ, 0.5; check_args=check_args)
# end

# DDM(ν::Real, α::Real, τ::Real, z::Real; check_args::Bool=true) = DDM(promote(ν, α, τ, z)...; check_args=check_args)
# DDM(ν::Real, α::Real, τ::Real; check_args::Bool=true) = DDM(promote(ν, α, τ)...; check_args=check_args)

I don't think you will be able to use DDM with a keyword argument because Julia unfortunately does not dispatch keyword arguments. A different approach will be needed if you would like to include checks on the correctness of the arguments. It might not be possible to opt in/out.

I'm not much of a fan about the negative RT scheme. It might be a clever way to make the model faster but I think it has two drawbacks worth considering. First, it has the potential to lead to confusion or errors. What if a person computes a mean over rts, or some other operation unaware of the presence of negative rts? Perhaps it might be easier to code correctly this way, but we only have to do that once. By contrast, users will have to use it appropriately everytime. A second drawback is that it is an exception to the general interface, which is to return an index vector of choices, and a vector of RTs. I think having a standardized interface is important because it sets up expectations for how the package works, and allows the user to potentially use the same code without having to modify depending on whether it uses LBA or DDM.

Another issue I noticed was that pdf returns NaN for negative RTs. I'm not sure why it does that. As a user, I would expect that I could pass simulated data to pdf without a problem. A similar problem exists for logpdf

Example

Random.seed!(5484)
dist = DDM(.5, .8, .3, .4)

rts = rand(dist, 10)

logpdf.(dist, rts)

Output

10-element Vector{Float64}:
 NaN
 NaN
  -2.0485034773320545
 NaN
 NaN
 NaN
  -2.6704531377602465
  -3.9876172218519876
 NaN
   0.34653192776453606

My guess is that there is a bug due to method ambiguity:

pdf(d::DDM{T}, t::Real; ϵ::Real = 1.0e-12) where {T<:Real}
pdf(d::DDM, t::Real; ϵ::Real = 1.0e-12)

I don't believe these methods can be distinguished based on their signatures, which causes the last one specified to be used. The last method has

    if τ ≥ t
        return T(NaN)
    end

which presumably is producing the NaN above.

But those are small details which should be easy to fix.

I am opposed to breaking the interface, but I wonder whether a solution would be to transform the data coming in and out? For example, in the pdf function, you might have:

rts = copy(_rts_in)
rts[choice .== 1] .= -rts[choice .==1]

Its less than ideal because it is copying data, but I would be amenable to that solution because it preserves the interface. If that is not a viable solution for you, I think the best approach moving forward would be keep the packages separate. In my opinion, the benefit of having models within the same package is that they follow the same interface. Eventually, I or someone else will add the full DDM (where the negative RT scheme possibly does not apply).

itsdfish commented 1 year ago

Also, I want to apologize for bringing up this issue now. If I knew that approach outputs negative RTs, I would have raised the issue before you made the PR.

kiante-fernandez commented 1 year ago

All good. Thanks @itsdfish catching those things.

I suggested the negative RT solution because it is canonical. People who have used these models are aware of it and feel comfortable with it, which makes it less error-prone when implementing and requires less coding. It was widely agreed upon that the JAGS version used a more intuitive coding approach that was less error-prone. For instance in Turing code, if only positive numbers were used and there was no negative sign, when faced with lower bound responses, the user now needs to invert both the starting point and the drift rate. This doubles the code length and increases the probability of making an error.

Additionally, I thought it might aid in our challenges with the Tupled output issue. But all that said, I agree that having a coherent interface across all the models is important, and the negative RT solution does not work when you have more than two boundaries (although the negative solution works in all cases involving only two).

Lets stick with the standardized interface. I can make the changes throughout for you to review!

itsdfish commented 1 year ago

Thanks for your understanding. Sounds like a plan.

In case you were not following the thread on discourse, I found a similar issue with predict for the Wald and InverseGaussian (implemented in Distributions.jl). Neither of those distributions return a tuple, but instead follow the interface for ContinuousUnivariateDistribution as intended. Based on that issue, it's not clear to me whether the Tuple is the problem (or at least only problem).

itsdfish commented 1 year ago

I am closing this in anticipation of merging the PR. By the way, I misunderstood the nature of the problem with predict.