From 0b825f9239786a9d7d20526ee5d29509ec2aab06 Mon Sep 17 00:00:00 2001 From: itsdfish Date: Sat, 29 Jun 2024 19:49:07 -0400 Subject: [PATCH 1/6] simplify increment and simulate --- src/LCA.jl | 2 +- src/MDFT.jl | 2 +- src/RDM.jl | 10 ++++------ 3 files changed, 6 insertions(+), 8 deletions(-) diff --git a/src/LCA.jl b/src/LCA.jl index 24c89972..8cf9be5e 100644 --- a/src/LCA.jl +++ b/src/LCA.jl @@ -177,5 +177,5 @@ function simulate(model::AbstractLCA; Δt = 0.001, _...) push!(evidence, copy(x)) push!(time_steps, t) end - return time_steps, reduce(vcat, transpose.(evidence)) + return time_steps, stack(evidence, dims = 1) end diff --git a/src/MDFT.jl b/src/MDFT.jl index a7d197f0..44ca07aa 100644 --- a/src/MDFT.jl +++ b/src/MDFT.jl @@ -367,5 +367,5 @@ function simulate(model::MDFT, M::AbstractArray; Δt = 0.001, _...) push!(evidence, copy(x)) push!(time_steps, t) end - return time_steps, reduce(vcat, transpose.(evidence)) + return time_steps, stack(evidence, dims = 1) end diff --git a/src/RDM.jl b/src/RDM.jl index 645e830f..04cc58c5 100644 --- a/src/RDM.jl +++ b/src/RDM.jl @@ -187,20 +187,18 @@ function simulate(rng::AbstractRNG, model::AbstractRDM; Δt = 0.001) z = rand(rng, Uniform(0, A), n) α = k + A x = z - ϵ = fill(0.0, n) evidence = [deepcopy(x)] time_steps = [t] while all(x .< α) t += Δt - increment!(rng, model, x, ϵ, ν; Δt) + increment!(rng, model, x, ν; Δt) push!(evidence, deepcopy(x)) push!(time_steps, t) end - return time_steps, reduce(vcat, transpose.(evidence)) + return time_steps, stack(evidence, dims = 1) end -function increment!(rng::AbstractRNG, model::AbstractRDM, x, ϵ, ν; Δt) - ϵ .= rand(rng, Normal(0.0, 1.0), length(ν)) - x .+= ν * Δt + ϵ * √(Δt) +function increment!(rng::AbstractRNG, model::AbstractRDM, x, ν; Δt) + x .+= ν * Δt .+ rand(rng, Normal(0.0, √(Δt)), length(ν)) return nothing end From ee07d575e47cdfca9f1a07e9cbca34478bed43ec Mon Sep 17 00:00:00 2001 From: itsdfish Date: Sun, 30 Jun 2024 07:47:20 -0400 Subject: [PATCH 2/6] start MLBA --- src/LBA.jl | 16 +++++----- src/MLBA.jl | 55 +++++++++++++++++++++++++++++++++ src/SequentialSamplingModels.jl | 3 ++ src/type_system.jl | 8 +++++ 4 files changed, 74 insertions(+), 8 deletions(-) create mode 100644 src/MLBA.jl diff --git a/src/LBA.jl b/src/LBA.jl index a0f63b22..f12265ab 100644 --- a/src/LBA.jl +++ b/src/LBA.jl @@ -5,11 +5,11 @@ A model object for the linear ballistic accumulator. # Parameters -- `ν`: a vector of drift rates -- `σ`: a vector of drift rate standard deviation -- `A`: max start point -- `k`: A + k = b, where b is the decision threshold -- `τ`: a encoding-response offset +- `ν::Vector{T}`: a vector of drift rates +- `σ::Vector{T}`: a vector of drift rate standard deviation +- `A::T`: max start point +- `k::T`: A + k = b, where b is the decision threshold +- `τ::T`: an encoding-response offset # Constructors @@ -50,13 +50,13 @@ function LBA(ν, σ, A, k, τ) return LBA(ν, σ, A, k, τ) end +LBA(; τ = 0.3, A = 0.8, k = 0.5, ν = [2.0, 1.75], σ = fill(1.0, length(ν))) = + LBA(ν, σ, A, k, τ) + function params(d::LBA) return (d.ν, d.σ, d.A, d.k, d.τ) end -LBA(; τ = 0.3, A = 0.8, k = 0.5, ν = [2.0, 1.75], σ = fill(1.0, length(ν))) = - LBA(ν, σ, A, k, τ) - function select_winner(dt) if any(x -> x > 0, dt) mi, mv = 0, Inf diff --git a/src/MLBA.jl b/src/MLBA.jl new file mode 100644 index 00000000..fcdb6e11 --- /dev/null +++ b/src/MLBA.jl @@ -0,0 +1,55 @@ +""" + MLBA{T <: Real} <: AbstractMLBA + +# Fields + +- `ν::Vector{T}`: a vector of drift rates, which is a function of β₀, λₚ, λₙ, γ +- `β₀::T`: baseline input for drift rate +- `λₚ::T`: decay constant for attention weights of positive differences +- `λₙ::T`: decay constant for attention weights of negative differences +- `γ::T`: risk aversion exponent for subjective values +- `σ::Vector{T}`: a vector of drift rate standard deviation +- `A::T`: max start point +- `k::T`: A + k = b, where b is the decision threshold +- `τ::T`: an encoding-response offset + +# References + +Trueblood, J. S., Brown, S. D., & Heathcote, A. (2014). The multiattribute linear ballistic accumulator model of context effects in multialternative choice. Psychological Review, 121(2), 179. +""" +mutable struct MLBA{T <: Real} <: AbstractMLBA + ν::Vector{T} + β₀::T + λₚ::T + λₙ::T + γ::T + σ::Vector{T} + A::T + k::T + τ::T +end + +function MLBA(ν, β₀, λₚ, λₙ, γ, σ, A, k, τ) + _, β₀, λₚ, λₙ, γ, _, A, k, τ = promote(ν[1], β₀, λₚ, λₙ, γ, σ[1], A, k, τ) + ν = convert(Vector{typeof(k)}, ν) + σ = convert(Vector{typeof(k)}, σ) + return MLBA(ν, β₀, λₚ, λₙ, γ, σ, A, k, τ) +end + +MLBA(; + n_alternatives = 3, + ν = fill(0.0, n_alternatives), + β₀ = 1.0, + λₚ = 1.0, + λₙ = 1.0, + γ = 0.70, + τ = 0.3, + A = 0.8, + k = 0.5, + σ = fill(1.0, n_alternatives) +) = + MLBA(ν, β₀, λₚ, λₙ, γ, σ, A, k, τ) + +function params(d::AbstractMLBA) + return (d.ν, d.β₀, d.λₚ, d.λₙ, d.γ, d.σ, d.A, d.k, d.τ) +end diff --git a/src/SequentialSamplingModels.jl b/src/SequentialSamplingModels.jl index a4e53b1f..47032d9f 100644 --- a/src/SequentialSamplingModels.jl +++ b/src/SequentialSamplingModels.jl @@ -38,6 +38,7 @@ export AbstractCDDM export AbstractLBA export AbstractLCA export AbstractLNR +export AbstractMLBA export AbstractMDFT export AbstractPoissonRace export AbstractRDM @@ -53,6 +54,7 @@ export LBA export LCA export LNR export maaDDM +export MLBA export MDFT export PoissonRace export SSM1D @@ -103,4 +105,5 @@ include("poisson_race.jl") include("stDDM.jl") include("MDFT.jl") include("ClassicMDFT.jl") +include("MLBA.jl") end diff --git a/src/type_system.jl b/src/type_system.jl index ba06e964..d29198db 100644 --- a/src/type_system.jl +++ b/src/type_system.jl @@ -58,6 +58,13 @@ An abstract type for the lognormal race model """ abstract type AbstractLNR <: SSM2D end +""" + AbstractMLBA <: AbstractLBA + +An abstract type for the multi-attribute linear ballistic accumulator +""" +abstract type AbstractMLBA <: AbstractLBA end + """ AbstractLCA <: SSM2D @@ -157,6 +164,7 @@ function rand(rng::AbstractRNG, d::SSM2D, n_trials::Int; kwargs...) end return (; choice, rt) end + rand(d::SSM2D, n_trials::Int; kwargs...) = rand(Random.default_rng(), d, n_trials; kwargs...) From 647d0661e1534e5a6ca0f569ea16ba363d32f7f2 Mon Sep 17 00:00:00 2001 From: itsdfish Date: Sun, 30 Jun 2024 08:32:49 -0400 Subject: [PATCH 3/6] use 1.9 for benchmark_pr --- .github/workflows/benchmark_pr.yml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/.github/workflows/benchmark_pr.yml b/.github/workflows/benchmark_pr.yml index ce10458b..d8bbd4b6 100644 --- a/.github/workflows/benchmark_pr.yml +++ b/.github/workflows/benchmark_pr.yml @@ -16,7 +16,7 @@ jobs: - uses: actions/checkout@v2 - uses: julia-actions/setup-julia@v1 with: - version: "1.8" + version: "1.9" - uses: julia-actions/cache@v1 - name: Extract Package Name from Project.toml id: extract-package-name From 5155e3d56185219caaadf880af3c6a66fa75beeb Mon Sep 17 00:00:00 2001 From: itsdfish Date: Sun, 30 Jun 2024 18:56:56 -0400 Subject: [PATCH 4/6] add mlba functions and tests --- docs/src/mdft.md | 2 +- src/MLBA.jl | 93 +++++++++++++++++++++++++++++++++++++++++++++++ test/mlba.jl | 95 ++++++++++++++++++++++++++++++++++++++++++++++++ test/runtests.jl | 2 +- 4 files changed, 190 insertions(+), 2 deletions(-) create mode 100644 test/mlba.jl diff --git a/docs/src/mdft.md b/docs/src/mdft.md index 3a3fbc31..0c75dac9 100644 --- a/docs/src/mdft.md +++ b/docs/src/mdft.md @@ -154,7 +154,7 @@ M₂ = [ choices,rts = rand(dist, 10_000, M₂; Δt = .001) probs2 = map(c -> mean(choices .== c), 1:2) ``` -Here, we see that job `A` is prefered over job `B`. Also note, in the code block above, `rand` has a keyword argument `Δt` which controls the precision of the discrete approximation. The default value is `Δt = .001`. +Here, we see that job `A` is prefered over job `B`. Also note, in the code block above, `rand` has a keyword argument `Δt` which controls the precision of the time discrete approximation. The default value is `Δt = .001`. Next, we will simulate the choice between jobs `A`, `B`, and `S`. diff --git a/src/MLBA.jl b/src/MLBA.jl index fcdb6e11..295844a7 100644 --- a/src/MLBA.jl +++ b/src/MLBA.jl @@ -13,6 +13,23 @@ - `k::T`: A + k = b, where b is the decision threshold - `τ::T`: an encoding-response offset +# Constructors + + MLBA(ν, β₀, λₚ, λₙ, γ, σ, A, k, τ) + + MLBA(; + n_alternatives = 3, + ν = fill(0.0, n_alternatives), + β₀ = 1.0, + λₚ = 1.0, + λₙ = 1.0, + γ = 0.70, + τ = 0.3, + A = 0.8, + k = 0.5, + σ = fill(1.0, n_alternatives) + ) + # References Trueblood, J. S., Brown, S. D., & Heathcote, A. (2014). The multiattribute linear ballistic accumulator model of context effects in multialternative choice. Psychological Review, 121(2), 179. @@ -53,3 +70,79 @@ MLBA(; function params(d::AbstractMLBA) return (d.ν, d.β₀, d.λₚ, d.λₙ, d.γ, d.σ, d.A, d.k, d.τ) end + +rand(d::AbstractMLBA, M::AbstractArray) = rand(Random.default_rng(), d, M) + +function rand(rng::AbstractRNG, d::AbstractMLBA, M::AbstractArray) + compute_drift_rates!(d, M) + return rand(rng, d) +end + +rand(d::AbstractMLBA, n_trials::Int, M::AbstractArray) = rand(Random.default_rng(), d, n_trials, M) + +function rand(rng::AbstractRNG, d::AbstractMLBA, n_trials::Int, M::AbstractArray) + compute_drift_rates!(d, M) + return rand(rng, d, n_trials) +end + +function compute_drift_rates!(dist::AbstractMLBA, M::AbstractArray) + (; ν, β₀) = dist + n_options = length(ν) + ν .= β₀ + utilities = map(x -> compute_utility(dist, x), eachrow(M)) + for i ∈ 1:n_options + for j ∈ 1:n_options + i == j ? continue : nothing + ν[i] += compare(dist, utilities[i], utilities[j]) + end + end + return nothing +end + +function compute_weight(dist::AbstractMLBA, u1, u2) + (; λₙ, λₚ) = dist + λ = u1 ≥ u2 ? λₚ : λₙ + return exp(-λ * abs(u1 - u2)) +end + +function compute_utility(dist::AbstractMLBA, v) + (; γ) = dist + θ = atan(v[2], v[1]) + # x and y intercepts for line passing through v with slope -1 + a = sum(v) + u = fill(0.0, 2) + u[1] = a / (tan(θ)^γ + 1)^(1 / γ) + u[2] = a * (1 - (u[1] / a)^γ)^(1 / γ) + #u[2] = (a * tan(θ)) / (1 + tan(θ)^γ)^(1/ γ) + return u +end + +function compare(dist::AbstractMLBA, u1, u2) + v = 0.0 + for i ∈ 1:2 + v += compute_weight(dist, u1[i], u2[i]) * (u1[i] - u2[i]) + end + return v +end + +""" + simulate(rng::AbstractRNG, model::AbstractMLBA, M::AbstractArray; n_steps = 100) + +Returns a matrix containing evidence samples of the MLBA decision process. In the matrix, rows +represent samples of evidence per time step and columns represent different accumulators. + +# Arguments + +- `model::AbstractMLBA`: a subtype of AbstractMLBA + +# Keywords + +- `n_steps=100`: number of time steps at which evidence is recorded +""" +function simulate(rng::AbstractRNG, model::AbstractMLBA, M::AbstractArray; n_steps = 100) + compute_drift_rates!(model, M) + return simulate(rng, model; n_steps) +end + +simulate(model::AbstractMLBA, M::AbstractArray; n_steps = 100) = + simulate(Random.default_rng(), model, M; n_steps) diff --git a/test/mlba.jl b/test/mlba.jl new file mode 100644 index 00000000..62f37856 --- /dev/null +++ b/test/mlba.jl @@ -0,0 +1,95 @@ +@safetestset "MLBA Tests" begin + @safetestset "compute_utility" begin + @safetestset "1" begin + using SequentialSamplingModels + using SequentialSamplingModels: compute_utility + using Test + + model = MLBA(; γ = 1) + + @test compute_utility(model, [3, 1]) ≈ [3, 1] + @test compute_utility(model, [2, 2]) ≈ [2, 2] + end + + @safetestset "2" begin + using SequentialSamplingModels + using SequentialSamplingModels: compute_utility + using Test + + model = MLBA(; γ = 0.5) + + @test compute_utility(model, [3, 1]) ≈ [1.60770, 0.53590] atol = 1e-4 + @test compute_utility(model, [2, 2]) ≈ [1, 1] + end + + @safetestset "3" begin + using SequentialSamplingModels + using SequentialSamplingModels: compute_utility + using Test + + model = MLBA(; γ = 1.2) + + @test compute_utility(model, [3, 1]) ≈ [3.2828, 1.0943] atol = 1e-4 + @test compute_utility(model, [2, 2]) ≈ [2.2449, 2.2449] atol = 1e-4 + end + end + + @safetestset "compute_weight" begin + using SequentialSamplingModels + using SequentialSamplingModels: compute_weight + using Test + + model = MLBA(; λₚ = 1.5, λₙ = 2.5) + + @test compute_weight(model, 2, 3) ≈ 0.082085 atol = 1e-4 + @test compute_weight(model, 3, 2) ≈ 0.22313 atol = 1e-4 + end + + @safetestset "compare" begin + using SequentialSamplingModels + using SequentialSamplingModels: compare + using Test + + model = MLBA(; λₚ = 1.5, λₙ = 2.5) + + @test compare(model, [3, 1], [1, 3]) ≈ 0.086098 atol = 1e-4 + @test compare(model, [1, 3], [3, 1]) ≈ 0.086098 atol = 1e-4 + @test compare(model, [2, 4], [3, 5]) ≈ -0.16417 atol = 1e-4 + end + + @safetestset "compute_drift_rates" begin + @safetestset "1" begin + using SequentialSamplingModels + using SequentialSamplingModels: compute_drift_rates! + using Test + + model = MLBA(; λₚ = 1.5, λₙ = 2.5, β₀ = 1, γ = 1.2) + + stimuli = [ + 1 3 + 3 1 + 2 2 + ] + + compute_drift_rates!(model, stimuli) + @test model.ν ≈ [1.2269, 1.2269, 1.2546] atol = 1e-4 + end + + @safetestset "2" begin + using SequentialSamplingModels + using SequentialSamplingModels: compute_drift_rates! + using Test + + model = MLBA(; λₚ = 0.75, λₙ = 3, β₀ = 2, γ = 2) + + stimuli = [ + 1 1 + 4 2 + 2 3 + ] + + compute_drift_rates!(model, stimuli) + @test model.ν ≈ [1.9480, 3.0471, 3.3273] atol = 1e-4 + end + end +end diff --git a/test/runtests.jl b/test/runtests.jl index f28b8acb..e10c37e8 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -2,4 +2,4 @@ using SafeTestsets files = filter(f -> f ≠ "runtests.jl", readdir()) -include.(files) +include.(files) \ No newline at end of file From 1d51f2e5fcaf7b703c77c566904c5b997ac1fab7 Mon Sep 17 00:00:00 2001 From: dfish <9677184+itsdfish@users.noreply.github.com> Date: Sun, 30 Jun 2024 19:07:47 -0400 Subject: [PATCH 5/6] Apply suggestions from code review Co-authored-by: github-actions[bot] <41898282+github-actions[bot]@users.noreply.github.com> --- src/MLBA.jl | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/src/MLBA.jl b/src/MLBA.jl index 295844a7..aa7fac0f 100644 --- a/src/MLBA.jl +++ b/src/MLBA.jl @@ -78,7 +78,8 @@ function rand(rng::AbstractRNG, d::AbstractMLBA, M::AbstractArray) return rand(rng, d) end -rand(d::AbstractMLBA, n_trials::Int, M::AbstractArray) = rand(Random.default_rng(), d, n_trials, M) +rand(d::AbstractMLBA, n_trials::Int, M::AbstractArray) = + rand(Random.default_rng(), d, n_trials, M) function rand(rng::AbstractRNG, d::AbstractMLBA, n_trials::Int, M::AbstractArray) compute_drift_rates!(d, M) From f38aef4e5e8c57d99018ecb3db0c1f96f788332c Mon Sep 17 00:00:00 2001 From: itsdfish Date: Wed, 3 Jul 2024 09:18:46 -0400 Subject: [PATCH 6/6] finish mlba documentation and tests --- Project.toml | 2 +- README.md | 1 + docs/make.jl | 1 + docs/src/mlba.md | 215 +++++++++++++++++++++++++++++++++++++++++++++++ src/MLBA.jl | 126 ++++++++++++++++++++++++--- test/mlba.jl | 62 ++++++++++++++ test/runtests.jl | 2 +- 7 files changed, 395 insertions(+), 14 deletions(-) create mode 100644 docs/src/mlba.md diff --git a/Project.toml b/Project.toml index 9cdb99a5..91d343a9 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "SequentialSamplingModels" uuid = "0e71a2a6-2b30-4447-8742-d083a85e82d1" authors = ["itsdfish"] -version = "0.11.0" +version = "0.11.1" [deps] Distributions = "31c24e10-a181-5473-b8eb-7969acd0382f" diff --git a/README.md b/README.md index 57ed4329..375c7ebf 100644 --- a/README.md +++ b/README.md @@ -22,6 +22,7 @@ The following SSMs are supported : - [Log Normal Race](https://itsdfish.github.io/SequentialSamplingModels.jl/dev/lnr/) - [Multi-attribute Attentional Drift Diffusion](https://itsdfish.github.io/SequentialSamplingModels.jl/dev/maaDDM/) - [Multi-attribute Decision Field Theory](https://itsdfish.github.io/SequentialSamplingModels.jl/dev/mdft/) +- [Multi-attribute Linear Ballistic Accumulator](https://itsdfish.github.io/SequentialSamplingModels.jl/dev/mlba/) - [Poisson Race](https://itsdfish.github.io/SequentialSamplingModels.jl/dev/poisson_race) - [Racing Diffusion](https://itsdfish.github.io/SequentialSamplingModels.jl/dev/rdm/) - [Wald](https://itsdfish.github.io/SequentialSamplingModels.jl/dev/wald/) diff --git a/docs/make.jl b/docs/make.jl index 982c648c..71f4d49d 100644 --- a/docs/make.jl +++ b/docs/make.jl @@ -32,6 +32,7 @@ makedocs( "Lognormal Race Model (LNR)" => "lnr.md", "Muti-attribute Attentional Drift Diffusion Model" => "maaDDM.md", "Multi-attribute Decision Field Theory" => "mdft.md", + "Multi-attribute Linear Ballistic Accumulator" => "mlba.md", "Poisson Race" => "poisson_race.md", "Racing Diffusion Model (RDM)" => "rdm.md", "Starting-time Drift Diffusion Model (stDDM)" => "stDDM.md", diff --git a/docs/src/mlba.md b/docs/src/mlba.md new file mode 100644 index 00000000..f23d5ccd --- /dev/null +++ b/docs/src/mlba.md @@ -0,0 +1,215 @@ +```@setup MLBA +using SequentialSamplingModels +using Plots +using Random +M = [ + 1.0 3.0 # A + 3.0 1.0 # B + 0.9 3.1 # S +] +``` + +# Multi-attribute Linear Ballistic Accumulator + +Multi-attribute Linear Ballistic Accumulator (MLBA; Trueblood et al., 2014) is an extention of the [LBA](https://itsdfish.github.io/SequentialSamplingModels.jl/dev/lba/) for multi-attribute deicisions. Alternatives in multi-attribute decisions vary along multiple dimensions. For example, jobs may differ in terms of benefits, salary, flexibility, and work-life balance. As with other sequential sampling models, MLBA assumes that evidence (or preference) accumulates dynamically until the evidence for one alternative reaches a threshold, and triggers the selection of the winning alternative. MLBA incorporates three additional core assumptions: + +1. Drift rates based on a weighted sum of pairwise comparisons +2. The comparison weights are an inverse function of similarity of two alternatives on a given attribute. +3. Objective values are mapped to subjective values, which can display extremeness aversion + +As with [MDFT](https://itsdfish.github.io/SequentialSamplingModels.jl/dev/mdft/), the MLBA excells in accounting for context effects in preferential decision making. A context effect occurs when the preference relationship between two alternatives changes when a third alternative is included in the choice set. In such cases, the preferences may reverse or the decision maker may violate rational choice principles. + +# Similarity Effect + +In what follows, we will illustrate the use of MLBA with a demonstration of the similarity effect. +Consider the choice between two jobs, `A` and `B`. The main criteria for evaluating the two jobs are salary and flexibility. Job `A` is high on salary but low on flexibility, whereas job `B` is low on salary. In the plot below, jobs `A` and `B` are located on the line of indifference, $y = 3 - x$. However, because salary recieves more attention, job `A` is slightly prefered over job `B`. + +```@example MLBA +scatter( + M[:, 1], + M[:, 2], + grid = false, + leg = false, + lims = (0, 4), + xlabel = "Flexibility", + ylabel = "Salary", + markersize = 6, + markerstrokewidth = 2 +) +annotate!(M[1, 1] + 0.10, M[1, 2] + 0.25, "A") +annotate!(M[2, 1] + 0.10, M[2, 2] + 0.25, "B") +annotate!(M[3, 1] + 0.10, M[3, 2] + 0.25, "S") +plot!(0:0.1:4, 4:-0.1:0, color = :black, linestyle = :dash) +``` +Suppose an job `S`, which is similar to A is added to the set of alternatives. Job `S` inhibits job `A` more than job `B` because `S` and `A` are close in attribute space. As a result, the preference for job `A` over job `B` is reversed. Formally, this is stated as: + +```math +\Pr(A \mid \{A,B\}) > \Pr(B \mid \{A,B\}) +``` + +```math +\Pr(A \mid \{A,B,S\}) < \Pr(B \mid \{A,B,S\}) +``` + +## Load Packages +The first step is to load the required packages. + +```@example MLBA +using SequentialSamplingModels +using Plots +using Random + +Random.seed!(8741) +``` +## Create Model Object +In the code below, we will define parameters for the MLBA and create a model object to store the parameter values. + +### Drift Rate Parameters +In sequential sampling models, the drift rate is the average speed with which evidence accumulates towards a decision threshold. In the MLBA, the drift rate is determined by comparing attributes between alternatives. The drift rate for alternative $i$ is defined as: + +$\nu_i = \beta_0 + \sum_{i\ne j} v_{i,j},$ + +where $\beta_0$ is the basline additive constant, and $v_{i,j}$ is the comparative value between alternatives $i$ and $j$. A given alternative $i$ is compared to another alternative $j \ne i$ as a weighted sum of differences across attributes $k \in [1,2,\dots n_a]$: + +$v_{i,j} = \sum_{k=1}^{n_a} w_{i,j,k} (u_{i,k} - u_{j,k}).$ +The attention weight between alternatives $i$ and $j$ on attribute $k$ is an inverse function of similarity, and decays exponentially with distance: + +$w_{i,j,k} = e^{-\lambda |u_{i,k} - u_{j,k}|}$ + +Similarity between alternatives is not necessarily symmetrical, giving rise to two decay rates: +```math +\lambda = \begin{cases} + \lambda_p & u_{i,k} \geq u_{j,k}\\ + \lambda_n & \mathrm{ otherwise}\\ +\end{cases} +``` +The subjective value $\mathbf{u} = [u_1,u_2]$ is found by bending the line of indifference passing through the objective stimulus $\mathbf{s}$ in attribute space, such that $(\frac{x}{a})^\gamma + (\frac{x}{b})^\gamma = 1$. When $\gamma > 1$, the model produces extremeness aversion. +#### Baseline Drift Rate + +The baseline drift rate parameter $\beta_0$ is a constant added to drift rate: + +```@example MLBA +β₀ = 5.0 +``` + +#### Similarity-Based Weighting + +Attention weights are an inverse function of similarity between alternatives on a given attribute. The decay rate for positive differences is: + +```@example MLBA +λₚ = 0.20 +``` +The decay rate for negative differences is: +```@example MLBA +λₙ = 0.40 +``` +The comparisons are asymmetrical when $\lambda_p \ne \lambda_n$. + +### Extremeness Aversion +The parameter for extremeness aversion is: +```@example MLBA +γ = 5 +``` + + $\gamma = 1$ indicates objective treatment of stimuli, whereas $\gamma > 1$ indicates extreness aversion, i.e. $[2,2]$ is prefered to $[3,1]$ even though both fall along the line indifference. + +### Standard Deviation of Drift Rates +The standard deviation of the drift rate distribution is given by $\sigma$, which is commonly fixed to 1 for each accumulator. + +```@example MLBA +σ = [1.0,1.0] +``` + +### Maximum Starting Point + +The starting point of each accumulator is sampled uniformly between $[0,A]$. + +```@example MLBA +A = 1.0 +``` +### Threshold - Maximum Starting Point + +Evidence accumulates until accumulator reaches a threshold $\alpha = k +A$. The threshold is parameterized this way to faciliate parameter estimation and to ensure that $A \le \alpha$. + +```@example MLBA +k = 1.0 +``` +### Non-Decision Time + +Non-decision time is an additive constant representing encoding and motor response time. + +```@example MLBA +τ = 0.30 +``` + +### MLBA Constructor + +Now that values have been asigned to the parameters, we will pass them to `MLBA` to generate the model object. We will begin with the choice between job `A` and job `B`. + +```@example MLBA +dist = MLBA(; + n_alternatives = 2, + β₀, + λₚ, + λₙ, + γ, + τ, + A, + k, +) +``` +## Simulate Model + +Now that the model is defined, we will generate 10,000 choices and reaction times using `rand`. + + ```@example MLBA +M₂ = [ + 1.0 3.0 # A + 3.0 1.0 # B +] + +choices,rts = rand(dist, 10_000, M₂) +probs2 = map(c -> mean(choices .== c), 1:2) +``` +Here, we see that job `A` is prefered over job `B`. + +Next, we will simulate the choice between jobs `A`, `B`, and `S`. + + ```@example MLBA +dist = MLBA(; + n_alternatives = 3, + β₀, + λₚ, + λₙ, + γ, + τ, + A, + k, +) + +M₃ = [ + 1.0 3.0 # A + 3.0 1.0 # B + 0.9 3.1 # S +] + +choices,rts = rand(dist, 10_000, M₃) +probs3 = map(c -> mean(choices .== c), 1:3) +``` +In this case, the preferences have reversed: job `B` is now preferred over job `A`. + +## Compute Choice Probability +The choice probability $\Pr(C=c)$ is computed by passing the model and choice index to `cdf` along with a large value for time as the second argument. + + ```@example MLBA +cdf(dist, 1, Inf, M₃) +``` + +## Plot Simulation +The code below plots a histogram for each alternative. + ```@example MLBA +histogram(dist; model_args = (M₃,)) +``` +# References + +Trueblood, J. S., Brown, S. D., & Heathcote, A. (2014). The multiattribute linear ballistic accumulator model of context effects in multialternative choice. Psychological Review, 121(2), 179. diff --git a/src/MLBA.jl b/src/MLBA.jl index aa7fac0f..af60ff80 100644 --- a/src/MLBA.jl +++ b/src/MLBA.jl @@ -20,13 +20,13 @@ MLBA(; n_alternatives = 3, ν = fill(0.0, n_alternatives), - β₀ = 1.0, - λₚ = 1.0, - λₙ = 1.0, - γ = 0.70, + β₀ = 5.0, + λₚ = 0.20, + λₙ = 0.40, + γ = 5.0, τ = 0.3, - A = 0.8, - k = 0.5, + A = 1.0, + k = 1.0, σ = fill(1.0, n_alternatives) ) @@ -56,13 +56,13 @@ end MLBA(; n_alternatives = 3, ν = fill(0.0, n_alternatives), - β₀ = 1.0, - λₚ = 1.0, - λₙ = 1.0, - γ = 0.70, + β₀ = 5.0, + λₚ = 0.20, + λₙ = 0.40, + γ = 5.0, τ = 0.3, - A = 0.8, - k = 0.5, + A = 1.0, + k = 1.0, σ = fill(1.0, n_alternatives) ) = MLBA(ν, β₀, λₚ, λₙ, γ, σ, A, k, τ) @@ -73,6 +73,16 @@ end rand(d::AbstractMLBA, M::AbstractArray) = rand(Random.default_rng(), d, M) +""" + rand(rng::AbstractRNG, d::AbstractMLBA, M::AbstractArray) + +Generates a single choice-rt pair of simulated data from the Multi-attribute Linear Ballistic Accumulator. + +# Arguments + +- `dist::AbstractMLBA`: an object for the multi-attribute linear ballistic accumulator +- `M::AbstractArray`: an alternative × attribute value matrix representing the value of the stimuli +""" function rand(rng::AbstractRNG, d::AbstractMLBA, M::AbstractArray) compute_drift_rates!(d, M) return rand(rng, d) @@ -81,11 +91,31 @@ end rand(d::AbstractMLBA, n_trials::Int, M::AbstractArray) = rand(Random.default_rng(), d, n_trials, M) +""" + rand(rng::AbstractRNG, d::AbstractMLBA, n_trials::Int, M::AbstractArray) + +Generates `n_trials` choice-rt pair of simulated data from the Multi-attribute Linear Ballistic Accumulator. + +# Arguments + +- `dist::AbstractMLBA`: an object for the multi-attribute linear ballistic accumulator +- `n_trials::Int`: the number of trials to simulate +- `M::AbstractArray`: an alternative × attribute value matrix representing the value of the stimuli +""" function rand(rng::AbstractRNG, d::AbstractMLBA, n_trials::Int, M::AbstractArray) compute_drift_rates!(d, M) return rand(rng, d, n_trials) end +""" + compute_drift_rates!(dist::AbstractMLBA, M::AbstractArray) + + +# Arguments + +- `dist::AbstractMLBA`: an object for the multi-attribute linear ballistic accumulator +- `M::AbstractArray`: an alternative × attribute value matrix representing the value of the stimuli +""" function compute_drift_rates!(dist::AbstractMLBA, M::AbstractArray) (; ν, β₀) = dist n_options = length(ν) @@ -100,12 +130,35 @@ function compute_drift_rates!(dist::AbstractMLBA, M::AbstractArray) return nothing end +""" + compute_weight(dist::AbstractMLBA, u1, u2) + +Computes attention weight between two options. The weight is an inverse function of similarity between the +options. Similarity between two options decays exponentially as a function of distance and is asymmetrical +when decay rates λₚ ≠ λₙ. + +# Arguments + +- `dist::AbstractMLBA`: an object for the multi-attribute linear ballistic accumulator +- `u1`: utility of the first option +- `u1`: utility of the second option +""" function compute_weight(dist::AbstractMLBA, u1, u2) (; λₙ, λₚ) = dist λ = u1 ≥ u2 ? λₚ : λₙ return exp(-λ * abs(u1 - u2)) end +""" + compute_utility(dist::AbstractMLBA, v) + +Computes the utility of an alternative. + +# Arguments + +- `dist::AbstractMLBA`: an object for the multi-attribute linear ballistic accumulator +- `v`: a vector representing the attributes of a given n_alternative +""" function compute_utility(dist::AbstractMLBA, v) (; γ) = dist θ = atan(v[2], v[1]) @@ -118,6 +171,18 @@ function compute_utility(dist::AbstractMLBA, v) return u end +""" + compare(dist::AbstractMLBA, u1, u2) + +Computes a weighted difference between two options, where the weight is an inverse function of similarity +between options. + +# Arguments + +- `dist::AbstractMLBA`: an object for the multi-attribute linear ballistic accumulator +- `u1`: utility of the first option +- `u1`: utility of the second option +""" function compare(dist::AbstractMLBA, u1, u2) v = 0.0 for i ∈ 1:2 @@ -126,6 +191,43 @@ function compare(dist::AbstractMLBA, u1, u2) return v end +""" + pdf(d::AbstractMLBA, c::Int, rt::Real, M::AbstractArray) + +Computes default probability density for multi-alternative linear ballistic accumulator. + +# Arguments + +- `dist::AbstractMLBA`: an object for the multi-attribute linear ballistic accumulator +- `c::Int`: choice index +- `rt::Real`: reaction time in seconds +- `M::AbstractArray`: an alternative × attribute value matrix representing the value of the stimuli +""" +function pdf(d::AbstractMLBA, c::Int, rt::Real, M::AbstractArray) + compute_drift_rates!(d, M) + return pdf(d, c, rt) +end + +""" + logpdf(d::AbstractMLBA, c::Int, rt::Real, M::AbstractArray) + +Computes default log probability density for multi-alternative linear ballistic accumulator. + +# Arguments + +- `dist::AbstractMLBA`: an object for the multi-attribute linear ballistic accumulator +- `c::Int`: choice index +- `rt::Real`: reaction time in seconds +- `M::AbstractArray`: an alternative × attribute value matrix representing the value of the stimuli +""" +function logpdf(d::AbstractMLBA, c::Int, rt::Real, M::AbstractArray) + compute_drift_rates!(d, M) + return logpdf(d, c, rt) +end + +loglikelihood(d::AbstractMLBA, data::NamedTuple, M::AbstractArray) = + sum(logpdf.(d, data..., (M,))) + """ simulate(rng::AbstractRNG, model::AbstractMLBA, M::AbstractArray; n_steps = 100) diff --git a/test/mlba.jl b/test/mlba.jl index 62f37856..d2c3dbea 100644 --- a/test/mlba.jl +++ b/test/mlba.jl @@ -92,4 +92,66 @@ @test model.ν ≈ [1.9480, 3.0471, 3.3273] atol = 1e-4 end end + + @safetestset "logpdf" begin + using Random + using SequentialSamplingModels + using SequentialSamplingModels: compute_drift_rates! + using Test + + Random.seed!(80071) + + parms = ( + λₚ = 0.20, + λₙ = 0.40, + β₀ = 5, + γ = 5, + τ = 0.3, + A = 0.8, + k = 0.5 + ) + + M = [ + 1 4 + 2 2 + 4 1 + ] + + choice, rt = rand(MLBA(; parms...,), 10000, M) + + λₚs = range(0.80 * parms.λₚ, 1.2 * parms.λₚ, length = 100) + LLs = map(λₚ -> sum(logpdf.(MLBA(; parms..., λₚ), choice, rt, (M,))), λₚs) + _, mxidx = findmax(LLs) + @test λₚs[mxidx] ≈ parms.λₚ rtol = 0.02 + + λₙs = range(0.80 * parms.λₙ, 1.2 * parms.λₙ, length = 100) + LLs = map(λₙ -> sum(logpdf.(MLBA(; parms..., λₙ), choice, rt, (M,))), λₙs) + _, mxidx = findmax(LLs) + @test λₙs[mxidx] ≈ parms.λₙ rtol = 0.02 + + β₀s = range(0.80 * parms.β₀, 1.2 * parms.β₀, length = 100) + LLs = map(β₀ -> sum(logpdf.(MLBA(; parms..., β₀), choice, rt, (M,))), β₀s) + _, mxidx = findmax(LLs) + @test β₀s[mxidx] ≈ parms.β₀ rtol = 0.02 + + γs = range(0.80 * parms.γ, 1.2 * parms.γ, length = 100) + LLs = map(γ -> sum(logpdf.(MLBA(; parms..., γ), choice, rt, (M,))), γs) + _, mxidx = findmax(LLs) + @test γs[mxidx] ≈ parms.γ rtol = 5e-2 + + τs = range(0.80 * parms.τ, 1.2 * parms.τ, length = 100) + LLs = map(τ -> sum(logpdf.(MLBA(; parms..., τ), choice, rt, (M,))), τs) + _, mxidx = findmax(LLs) + @test τs[mxidx] ≈ parms.τ rtol = 0.02 + + As = range(0.80 * parms.A, 1.2 * parms.A, length = 100) + LLs = map(A -> sum(logpdf.(MLBA(; parms..., A), choice, rt, (M,))), As) + _, mxidx = findmax(LLs) + @test As[mxidx] ≈ parms.A rtol = 0.02 + + ks = range(0.80 * parms.k, 1.2 * parms.k, length = 100) + LLs = map(k -> sum(logpdf.(MLBA(; parms..., k), choice, rt, (M,))), ks) + _, mxidx = findmax(LLs) + @test ks[mxidx] ≈ parms.k rtol = 0.02 + end end diff --git a/test/runtests.jl b/test/runtests.jl index e10c37e8..f28b8acb 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -2,4 +2,4 @@ using SafeTestsets files = filter(f -> f ≠ "runtests.jl", readdir()) -include.(files) \ No newline at end of file +include.(files)