Should vectorization affect the adaptation or sampling paths given the same seed?

General
1

I have a model where I do something like (dumbed down):

parameters{
  vector[nrow_y] mu ;
}
model{
  for(col in 1:ncol_y){
    for(row in 1:nrow_y){
      y[row,col] ~ normal(mu[row],1) ;
    }
  }
}

Which of course can be sped up by vectorizing:

parameters{
  vector[n] mu ;
}
model{
  for(col in 1:ncol_y){
    y[,col] ~ normal(mu,1) ;
  }
)

Using cmdstanr and keeping the data and seed constant, I thought that the only consequence of switching from the first to the second model would be to speed up computation, but that the actual course of adaptation and sampling would remain the same since the computations should yield the same values. However, this seems not to be the case as running with the same data and seed constant I get entirely different values in the posterior (and noticeably different adaptation performance).

Now, this arose while debugging a more complex model that definitely has it’s own pathologies I’m still working out, so maybe I’ve done something else wrong leading to the non-equivalence, but I wanted to check first: am I wrong to even expect equivalence? Does the vectorization change something about the compute state that would lead to divergence of the adaptation paths?

2

In theory, no, but in practice possibly due to floating point reasons.

2 Likes