Different outcomes with the same function

General Usage Performance
1

Hi, I’m making a program in Julia that calculates the posteriori distribution of a beta binomial model.

My first model was the following:

@model betabinomial(n, z, x, y, delta, gamma) = begin
alpha_p ~ truncated(Normal(x, 0.1), 0, Inf)
beta_p ~ truncated(Normal(y, 0.1), 0, Inf)
q_ajustada = Vector{Real}(undef, 3)
p = Vector{Real}(undef, 3)
for i = 1:length(n)
p[i] ~ Beta(alpha_p, beta_p)
q_ajustada[i] = deltap[i] + (1-gamma)(1-p[i])
z[i] ~ Binomial(n[i], q_ajustada[i])
end
end;

But now, instead of using the for loop, I want to use vector to achieve a better performance and less computing time so my new code looks like this:

@model betabinomial_v2(n::Array{Int64, 1}, z::Array{Int64, 1}, x::Float64, y::Float64, delta::Float64, gamma::Float64) = begin
#Definimos las iniciales
alpha_p ~ truncated(Normal(x, 0.1), 0, Inf)
beta_p ~ truncated(Normal(y, 0.1), 0, Inf)

# p se distribuye Beta
p ~ filldist(Beta(alpha_p, beta_p), length(n))

#q_ajustada sí tiene que ser así segun yo
q_ajustada = Vector{Float64}(undef, length(n))
unos = ones(length(n))
q_ajustada = delta*p + (1 - gamma)*(unos - p)

@. z ~ Binomial(n, q_ajustada)

end;

I defined the type of each variable and follow this publication Vectorising Turing models for the vectorization.
However, even when I’m using a seed, I don’t get the same results using both models. Any ideas of why and if I’m making a mistake, please let me know.

Thank you in advanced :slight_smile:

2

First of all, welcome!

Second, before you get too far into thinking about the new version of your code, let’s start with your premise:

This simply isn’t true in Julia. There’s no reason to prefer “vectorized” code in this form for speed compared to code using loops. Loops in Julia are as fast or faster than vectorized code in slow languages like python or matlab. If you’re happy with your loop-based code, then you’re done :slightly_smiling_face:

That said, there are some easy things you can do to improve performance of your original code. You have:

which is a container (Vector) with an abstract type (Real). This is exactly what’s warned against here: Performance Tips · The Julia Language

You can probably use something like Vector{typeof(alpha_p)} to improve performance without needing to rewrite the algorithm. Use @code_warntype and check out the other performance tips from that document to learn more about writing faster code.

2 Likes