Initialization of unitvector in rstan - confusing error message

Hi,

i tried to include init parameters to my model via the rstan interface. I have a unitvector unit_vector[2] v_a_mu[K] and filled it with values of the sine and cosine of K angles.

I defined my init function like this:

init_fun <- function(){
v_a_mu <- array(NA, dim = c(K, 2))
kapp <- NULL
d_mu <- NULL
d_sigma <- NULL

values <- list(v_a_mu = v_a_mu, d_mu = d_mu, kapp = kapp, d_sigma = d_sigma, rho_1 = rho_1, rho_2 = rho_2) 

return(list(values))
}

and defined my arrays, including v_a_mu, inside the function like this:

v_a_mu[1,1] <- -0.416146836547
v_a_mu[1,2] <- 0.909297426826

v_a_mu[2,1] <- 0.540302305868
v_a_mu[2,2] <- 0.841470984808

v_a_mu[3,1] <- -0.540302305868
v_a_mu[3,2] <- -0.841470984808

v_a_mu[4,1] <- -0.9899924966 
v_a_mu[4,2] <- 0.14112000806

v_a_mu[5,1] <- -0.9899924966
v_a_mu[5,2] <- -0.14112000806 

v_a_mu[6,1] <- -0.9899924966
v_a_mu[6,2] <- -0.14112000806

v_a_mu[7,1] <- -0.416146836547
v_a_mu[7,2] <- 0.909297426826

This didn’t work, and didn’t work as well when I defined init_r to [-pi, pi], the range of my angles.

Rejecting initial value:
Log probability evaluates to log(0), i.e. negative infinity.
Stan can't start sampling from this initial value.
Initialization between (-3.14159, 3.14159) failed after 100 attempts. 
 Try specifying initial values, reducing ranges of constrained values, or reparameterizing the model.

My model works fine if I don"t specify an init, so I thought I’d take the values of v_a_mu which are produced by my stan model when I don’t define a specific init, because with those sampling can be done. I printed them and filled my unitvector

v_a_mu[1,1] <- 0.978121
v_a_mu[1,2] <- 0.208038

v_a_mu[2,1] <- -0.742142
v_a_mu[2,2] <- 0.670242

v_a_mu[3,1] <- 0.624176
v_a_mu[3,2] <- 0.781283

v_a_mu[4,1] <- 0.850055 
v_a_mu[4,2] <- 0.526694

v_a_mu[5,1] <- -0.254035
v_a_mu[5,2] <- -0.967195 

v_a_mu[6,1] <- -0.42736
v_a_mu[6,2] <- -0.904081

v_a_mu[7,1] <- 0.761295
v_a_mu[7,2] <- -0.648406

And I call my model with

uni <- stan(uni_mu.stan", data = list(K = K, N = N, angle = angle, dist = dist, d_min = d_min), chains = nb_chain, iter = iteration, init = init_fun(), init_r = pi ,cores = 4)

but I get this error messages

SAMPLING FOR MODEL 'uni_mu' NOW (CHAIN 1).
[1] "Error in sampler$call_sampler(args_list[[i]]) : "

 [2] "  Error transforming variable v_a_mu: stan::io::unit_vector_unconstrain: Vector is not a valid unit vector. The sum of the squares of the elements should be 1, but is 1"
error occurred during calling the sampler; sampling not done
Stan model 'uni_mu' does not contain samples.

I don’t know what my mistake is, but this second error message doesn’t make any sense?

I guess something has to be wrong with my list definition, i can’t think of anything else. Does anyone have an idea?

Hi,

Some comments below as to what might be going wrong:

  1. K and rho_* are not defined in your init_fun.
  2. The parameters are not actually being initialized in init_fun (they are NULL).
  3. init_fun returns a list of a list when it should just return a list.

The rstan doc on the stan function (?rstan::stan) has an example that you can look at.

edit: Not a problem! I messed up, see post below!

Just checked, I’m pretty sure this is a problem. I went ahead a made an issue for it: https://github.com/stan-dev/stan/issues/2366. Thanks for bringing this up Jennifer!

Regarding your model, if Stan is able to initialize this model without specifying initial conditions, you should honestly probly go with that (so you don’t have to worry about your initial conditions biasing your results).

Wait, no, I messed up. The vectors you’re passing in aren’t unit length by floating point standards. It’s an issue with text representation of real numbers vs. their representation as doubles in the computer.

For this data:

v_a_mu[7,1] <- 0.761295
v_a_mu[7,2] <- -0.648406
sprintf("%.10f", 1.0 - sum(c(0.761295, -0.648406)^2))

Returns

"-0.0000004179"

If you’re happy with your vectors, set them and renormalize em’ before you go to use them.

Something like:

for(i in 1:7) {
  v_a_mu[i,] = v_a_mu[i, ] / sqrt(sum(va_a_mu[i,]^2))
}

(this even happens for the really high precision numbers you had in your first initialization):

sprintf("%.17f", 1.0 - sum(c(-0.416146836547, 0.909297426826)^2))

returns

"-0.00000000000046052"

edit: double precision goes out to like 17 digits or so I think

1 Like

Ah I see, haven’t thought about this… The error message still is very confusing :D but thank you very much!

There’s a discussion in the manual on the vagaries of floating point, but it’s only a sketch.

You should be getting an error message that the unit vector isn’t a unit vector. If you’re not seeing it, you may be using 2.14 or 2.15, where there was a bug suppressing these messages getting printed.

The initialization range is on the unconstrained scale. So if there are constraints, you’ll be transforming the unconstrained inits back to the constrained scales. For example, scaled and shifted inverse logit is used for variables declared with upper and lower bounds. If your parameters are bounded, you must use the appropriate constraints.

@bbbales2: Stan stops short of trying to validate constraints all the way down to floating point, but it’s pretty strict. I think most of the constraints are evaluated at around half the log precision, about 1e-8 or so.