TuringLang / ParetoSmooth.jl

An implementation of PSIS algorithms in Julia.
http://turinglang.org/ParetoSmooth.jl/
MIT License
19 stars 12 forks source link

loo_compare: UndefVarError: `Not` not defined #91

Closed smith-garrett closed 4 months ago

smith-garrett commented 4 months ago

I am trying to use the loo_compare function to compare two models in Turing, and I'm running into the error in the title. Here is a minimal working example (based on the official docs) where I run into the issue. (I'm running this in a Pluto notebook.):

using Distributions, Turing, ParetoSmooth

@model function model1(data)
        μ ~ Normal()
        for i in 1:length(data)
            data[i] ~ Normal(μ, 1)
        end
end

@model function model2(data)
        μ ~ Normal(1, 0.2)
        for i in 1:length(data)
            data[i] ~ Normal(μ, 1)
        end
end

data = rand(Normal(0, 1), 100)

chain = sample(model1(data), NUTS(), 1000)
chain2 = sample(model2(data), NUTS(), 1000)
p1 = psis_loo(model(data), chain)
p2 = psis_loo(model2(data), chain2)
loo_compare(p1, p2)

I have tried importing InvertedIndices explicitly, which seems to be where Not() comes from, but the error still happens. Any ideas on what might be going on? Thanks in advance!

Version info:

devmotion commented 4 months ago

I assume the problem is Not in https://github.com/TuringLang/ParetoSmooth.jl/blob/8a277cb45139510a51366f931ba9c4e8c86d13c1/src/ModelComparison.jl#L101. Unclear to me why one would use Not there instead of e.g. 2:end (if 1-based indexing is guaranteed which I did not check).

smith-garrett commented 4 months ago

Yes, that is the line that was referenced in the stack trace.

Glancing at the code, it looks like it should be safe to replace the Nots with 2:end, like you suggested. I'm happy to work on a pull request next week if it's actually that simple.

ParadaCarleton commented 4 months ago

Unclear to me why one would use Not there instead of e.g. 2:end (if 1-based indexing is guaranteed which I did not check).

Easier to read.

I don't know if it's guaranteed, but it's best not to tempt fate regardless. This is code from 4 years ago when I was still new to Julia and thought everything would always be 1-indexed.

smith-garrett commented 4 months ago

Yes, it's probably best to keep compatibility with different indexing schemes. I've never worked with OffsetArrays or the like, but one option I can see would be to index like this:

lastidx = axes(pointwise_diffs, 3)[2:end]
@views @. pointwise_diffs[:, 1:3, lastidx]...

This isn't as compact as Not, but it seems like we would need to introduce InvertedIndices as a dependency otherwise to make it work. Do you have a preference for either of the two options?

ParadaCarleton commented 4 months ago

InvertedIndices.jl is lightweight and going to be a dependency for basically any real data science workflow; feel free to make a PR!

devmotion commented 4 months ago

Fixed by #92